Tek-Tips is the largest IT community on the Internet today!

Members share and learn making Tek-Tips Forums the best source of peer-reviewed technical information on the Internet!

Register Log in
  • Congratulations strongm on being selected by the Tek-Tips community for having the most helpful posts in the forums last week. Way to Go!

subroutine errors

Status
Not open for further replies.
Tony1984

Tony1984

Technical User
Jul 2, 2012
20
US
Dear All
i found a subroutine in the google search for fitting nonliear data, but it seems the subroutine haa a lot of synatix errors, im not a proffessional in fortran
first: any one can help me to solve this errors ,
second: how any one help me to connect this subroutine to a program in order to make an .exe file
Thanks

SUBROUTINE LLSQFVL(NDATA,NPARM,MXPARM,YO,YU,YD,PV,PU,PS,CM,DSE)
INTEGER I,J,K,L,IDF,NDATA,NPARM,MXPARM
REAL*8 YO(NDATA), YU(NDATA), YD(NDATA), PV(NPARM), PU(NPARM),
PS(NPARM), CM(MXPARM,MXPARM), DSE,
PX(301), F95(10), TFACT, S, U
DATA F95/12.7062D0,4.3027D0,3.1824D0,2.7764D0,2.5706D0,2.4469D0,
2.3646D0,2.3060D0,2.2622D0,2.2281D0/
IF((NPARM.GT.MXPARM).OR.(NPARM.GT.301)
.OR.(NPARM.GT.NDATA)) THEN
WRITE(6,601) NDATA,NPARM,MXPARM
STOP
ENDIF
IF(NPARM.LE.0) THEN
DSE= 0.D0
DO 2 I= 1,NDATA
YD(I)= YO(I)
DSE= DSE+ (YD(I)/YU(I))**2
DSE= DSQRT(DSE/DFLOAT(NDATA))
RETURN
ENDIF
TFACT= 0.D0
IF(NDATA.GT.NPARM) THEN
IDF= NDATA-NPARM
IF(IDF.GT.10) TFACT= 1.960D0*DEXP(1.265D0/DFLOAT(IDF))
IF(IDF.LE.10) TFACT= F95(IDF)
ELSE
TFACT= 0.D0
ENDIF
DO 10 I = 1,NPARM
PS(I) = 0.D0
PU(I) = 0.D0
PX(I) = 0.D0
DO 8 J = 1,NPARM
CM(I,J) = 0.D0
CONTINUE
DO 14 I = 1,NDATA
S = 1.D0 / YU(I)
U = YO(I) * S
CALL DYIDPJ(I,PV)
DO 12 J = 1,NPARM
PV(J) = PV(J)*S
CALL QROD(NPARM,MXPARM,MXPARM,CM,PV,PX,U,PS,PU)
CONTINUE
CM(1,1) = 1.D0 / CM(1,1)
DO 20 I = 2,NPARM
L = I - 1
DO 18 J = 1,L
S = 0.D0
DO 16 K = J,L
S = S + CM(K,I) * CM(J,K)
CM(J,I) = -S / CM(I,I)
CM(I,I) = 1.D0 / CM(I,I)
DO 26 I = 1,NPARM
J = NPARM - I + 1
PV(J) = 0.D0
DO 24 K = J,NPARM
PV(J) = PV(J) + CM(J,K) * PX(K)
CONTINUE
DO 30 I = 1,NPARM
DO 30 J = I,NPARM
U = 0.D0
DO 28 K = J,NPARM
U = U + CM(I,K) * CM(J,K)
CM(I,J) = U
DO 36 J = 1,NPARM
PU(J) = DSQRT(CM(J,J))
DO 32 K= J,NPARM
CM(J,K)= CM(J,K)/PU(J)
DO 34 K= 1,J
CM(K,J)= CM(K,J)/PU(J)
CM(J,K)= CM(K,J)
PX(J)= 0.d0
DSE= 0.D0
DO 40 I = 1,NDATA
S = 1.D0 / YU(I)
U = 0.D0
CALL DYIDPJ(I,PS)
DO 38 J = 1,NPARM
PX(J)= PX(J) + (PS(J)*S)**2
U = U + PS(J) * PV(J)
YD(I) = YO(I) - U
DSE= DSE+ (S*YD(I))**2
IF(NDATA.GT.NPARM) THEN
DSE= DSQRT(DSE/(NDATA-NPARM))
ELSE
DSE= 0.d0
ENDIF
U= DSE*0.1d0/DFLOAT(NPARM)
S= DSE*TFACT
DO 44 J = 1,NPARM
PU(J)= S* PU(J)
PS(J)= U*DSQRT(NDATA/PX(J))
RETURN
601 FORMAT(/' *** Dimensioning problems in LLSQF *** (NDATA, NPARM, MX &
1PARM) = (',I5,2(' ,',I5),' )')
END

SUBROUTINE QROD(N,NR,NC,A,R,F,B,GC,GS)
INTEGER I,J,K,N,NC,NR
REAL*8 A(NR,NC), R(N), F(NR), GC(N), GS(N), B, Z(2)
DO 10 I = 1,N
Z(1) = R(I)
J = I - 1
DO 2 K = 1,J
Z(2) = GC(K) * A(K,I) + GS(K) * Z(1)
Z(1) = GC(K) * Z(1) - GS(K) * A(K,I)
A(K,I) = Z(2)
CONTINUE
GC(I) = 1.D0
GS(I) = 0.D0
IF(Z(1) .EQ. 0.D0) GOTO 10
IF(DABS(A(I,I)) .GE. DABS(Z(1))) GOTO 3
Z(2) = A(I,I) / Z(1)
GS(I) = 1.D0 / DSQRT(1.D0 + Z(2) * Z(2))
GC(I) = Z(2) * GS(I)
GOTO 4
Z(2) = Z(1) / A(I,I)
GC(I) = 1.D0 / DSQRT(1.D0 + Z(2) * Z(2))
GS(I) = Z(2) * GC(I)
A(I,I) = GC(I) * A(I,I) + GS(I) * Z(1)
Z(2) = GC(I) * F(I) + GS(I) * B
B = GC(I) * B - GS(I) * F(I)
F(I) = Z(2)
CONTINUE
RETURN
END

SUBROUTINE DYIDPJ(i,PD)
COMMON MXDATA,MXPARM,NPARM,X
REAL*8 PD(MXPARM),X(MXDATA),POWER
POWER= 1.D0
PD(1)= POWER
DO 10 J= 2,NPARM
POWER= POWER*X(I)
PD(J)= POWER
RETURN
END
 
Sort by date Sort by votes
Tony,

for this and your other thread it would be a good idea to enclose your code in code-tags. This would keep the original format and it would be easy to copy and paste it into my compilerwindow.

As far as I can see all the labels seem to be missing, but might have got lost with the format.

What compiler are you using ?

Fitting curves by least square method is not so uncommon a problem and there are quite a few softwarepackages out there, so it might be the best way to forget about this apparently faulty one and look for another in that link I gave you in your older thread.

Norbert





The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.
 
Upvote 0 Downvote
Code: 
SUBROUTINE LLSQFVL(NDATA,NPARM,MXPARM,YO,YU,YD,PV,PU,PS,CM,DSE)
INTEGER I,J,K,L,IDF,NDATA,NPARM,MXPARM
REAL*8 YO(NDATA), YU(NDATA), YD(NDATA), PV(NPARM), PU(NPARM),
PS(NPARM), CM(MXPARM,MXPARM), DSE,
PX(301), F95(10), TFACT, S, U
DATA F95/12.7062D0,4.3027D0,3.1824D0,2.7764D0,2.5706D0,2.4469D0,
2.3646D0,2.3060D0,2.2622D0,2.2281D0/
IF((NPARM.GT.MXPARM).OR.(NPARM.GT.301)
.OR.(NPARM.GT.NDATA)) THEN
WRITE(6,601) NDATA,NPARM,MXPARM
STOP
ENDIF
IF(NPARM.LE.0) THEN
DSE= 0.D0
DO 2 I= 1,NDATA
YD(I)= YO(I)
DSE= DSE+ (YD(I)/YU(I))**2
DSE= DSQRT(DSE/DFLOAT(NDATA))
RETURN
ENDIF
TFACT= 0.D0
IF(NDATA.GT.NPARM) THEN
IDF= NDATA-NPARM
IF(IDF.GT.10) TFACT= 1.960D0*DEXP(1.265D0/DFLOAT(IDF))
IF(IDF.LE.10) TFACT= F95(IDF)
ELSE
TFACT= 0.D0
ENDIF
DO 10 I = 1,NPARM
PS(I) = 0.D0
PU(I) = 0.D0
PX(I) = 0.D0
DO 8 J = 1,NPARM
CM(I,J) = 0.D0
CONTINUE
DO 14 I = 1,NDATA
S = 1.D0 / YU(I)
U = YO(I) * S
CALL DYIDPJ(I,PV)
DO 12 J = 1,NPARM
PV(J) = PV(J)*S
CALL QROD(NPARM,MXPARM,MXPARM,CM,PV,PX,U,PS,PU)
CONTINUE
CM(1,1) = 1.D0 / CM(1,1)
DO 20 I = 2,NPARM
L = I - 1
DO 18 J = 1,L
S = 0.D0
DO 16 K = J,L
S = S + CM(K,I) * CM(J,K)
CM(J,I) = -S / CM(I,I)
CM(I,I) = 1.D0 / CM(I,I)
DO 26 I = 1,NPARM
J = NPARM - I + 1
PV(J) = 0.D0
DO 24 K = J,NPARM
PV(J) = PV(J) + CM(J,K) * PX(K)
CONTINUE
DO 30 I = 1,NPARM
DO 30 J = I,NPARM
U = 0.D0
DO 28 K = J,NPARM
U = U + CM(I,K) * CM(J,K)
CM(I,J) = U
DO 36 J = 1,NPARM
PU(J) = DSQRT(CM(J,J))
DO 32 K= J,NPARM
CM(J,K)= CM(J,K)/PU(J)
DO 34 K= 1,J
CM(K,J)= CM(K,J)/PU(J)
CM(J,K)= CM(K,J)
PX(J)= 0.d0
DSE= 0.D0
DO 40 I = 1,NDATA
S = 1.D0 / YU(I)
U = 0.D0
CALL DYIDPJ(I,PS)
DO 38 J = 1,NPARM
PX(J)= PX(J) + (PS(J)*S)**2
U = U + PS(J) * PV(J)
YD(I) = YO(I) - U
DSE= DSE+ (S*YD(I))**2
IF(NDATA.GT.NPARM) THEN
DSE= DSQRT(DSE/(NDATA-NPARM))
ELSE
DSE= 0.d0
ENDIF
U= DSE*0.1d0/DFLOAT(NPARM)
S= DSE*TFACT
DO 44 J = 1,NPARM
PU(J)= S* PU(J)
PS(J)= U*DSQRT(NDATA/PX(J))
RETURN
601 FORMAT(/' *** Dimensioning problems in LLSQF *** (NDATA, NPARM, MX &
1PARM) = (',I5,2(' ,',I5),' )')
END

SUBROUTINE QROD(N,NR,NC,A,R,F,B,GC,GS)
INTEGER I,J,K,N,NC,NR
REAL*8 A(NR,NC), R(N), F(NR), GC(N), GS(N), B, Z(2)
DO 10 I = 1,N
Z(1) = R(I)
J = I - 1
DO 2 K = 1,J
Z(2) = GC(K) * A(K,I) + GS(K) * Z(1)
Z(1) = GC(K) * Z(1) - GS(K) * A(K,I)
A(K,I) = Z(2)
CONTINUE
GC(I) = 1.D0
GS(I) = 0.D0
IF(Z(1) .EQ. 0.D0) GOTO 10
IF(DABS(A(I,I)) .GE. DABS(Z(1))) GOTO 3
Z(2) = A(I,I) / Z(1)
GS(I) = 1.D0 / DSQRT(1.D0 + Z(2) * Z(2))
GC(I) = Z(2) * GS(I)
GOTO 4
Z(2) = Z(1) / A(I,I)
GC(I) = 1.D0 / DSQRT(1.D0 + Z(2) * Z(2))
GS(I) = Z(2) * GC(I)
A(I,I) = GC(I) * A(I,I) + GS(I) * Z(1)
Z(2) = GC(I) * F(I) + GS(I) * B
B = GC(I) * B - GS(I) * F(I)
F(I) = Z(2)
CONTINUE
RETURN
END

SUBROUTINE DYIDPJ(i,PD)
COMMON MXDATA,MXPARM,NPARM,X
REAL*8 PD(MXPARM),X(MXDATA),POWER
POWER= 1.D0
PD(1)= POWER
DO 10 J= 2,NPARM
POWER= POWER*X(I)
PD(J)= POWER
RETURN
END
 
Upvote 0 Downvote
I think there is something essential missing. This cannot be the code as you downloaded it.

Click the code con first, then include your code between the tags. The inputbox will modify the format of your code if you do it the other way round, that means reduces tabs to single blanks and I assume cut off the first column. You see, all the statement labels, that must have been present in the past, are not included.

A proper loop structure would be like

Code: 
    DO 10 J= 2,NPARM
    POWER= POWER*X(I)
10  PD(J)= POWER

to define start and end of the loop.
If your download looks as in your response, then forget this code. For you would have to add the labels to define loop ends and you do not know where, I assume.

Norbert




The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.
 
Upvote 0 Downvote
Hi Tony1984,
Nothing in wrong, but looking at your questions, I thing that you will not be productive with Fortran to solve your problem. It seems that you know nothing about Fortran programming and very little about fitting algorithms. I would suggest you to use a free system for scientific computing, where you don't have to bother with long program code.
For example:
* Matlab like systems: Scilab or Octave
* Statistical system: R
* Python for scientific computing SciPy, NumPy
etc ...

Here I found something about fitting in R, NumPy, Scilab and Octave

If you are interested, then google for more informations...
 
Upvote 0 Downvote
Status
Not open for further replies.

Similar threads

virk24
  • Locked
  • Question
subroutine call in a function
Replies
2
Views
306
virk24
virk24
drudox
  • Locked
  • Question
Subroutine give me very strange error
Replies
1
Views
97
xwb
xwb
FloatingFeather
  • Locked
  • Question
Fortran90 code with FORTRAN77 subroutine 1
Replies
3
Views
88
xwb
xwb
Mike Huffman
  • Locked
  • Question
core dump when returning from subroutine in fortran 77
Replies
4
Views
90
Mike Huffman
Mike Huffman
wassim2015
  • Locked
  • Question
subroutine does not run
Replies
1
Views
63
xwb
xwb

Part and Inventory Search

Sponsor

Back
Top