Comp Sci Fortran77 woes- single parameter least squares minimisation stuck

[NathanM]

Hey guys! Having issues with a code I'm writing. It doesn't fit the typical format and I may have already accidentally posted it in the wrong forum, so here's a copy of what was there with the newest version of the code. It'll (eventually) be a four-parameter least squares minimisation, but for the moment I'm just trying to get the single parameter version to work. I've got the first module done but, when I run it it gets caught up in a loop and doesn't actually minimise anything! Good start huh?

Here's my code, bless me with your insight!

PROGRAM assignment

! A program designed to fit experiemental data, using the method
! of least squares to minimise the associated chi-squared and
! obtain the four control parameters A,B,h1 and h2.
!*****************************************************************

IMPLICIT NONE
INTEGER i,t0
DOUBLE PRECISION t(17),t2(17),Ct(17),eCt(17),Ctdiff(17),c(17)
DOUBLE PRECISION Aa,Ab,Ba,Bb,h1a,h1b,h2a,h2b,chisqa,chisqb,dchisq
DOUBLE PRECISION deltah,Cs(17)

OPEN(21, FILE='data.txt', FORM='FORMATTED', STATUS='OLD')
DO i=1,17
READ(21,*)t(i),Ct(i),eCt(i)
END DO
CLOSE(21)

!Read in data.txt as three one dimensional arrays.
!*****************************************************************
!OPEN(21, FILE='outtest.txt', FORM='FORMATTED', STATUS='NEW')
!DO i=1,17
! WRITE(21,*)t(i),Ct(i),eCt(i)
!END DO
!CLOSE(21)
!
!Just to check input file is being read correctly.
!*****************************************************************
!****************************Parameters***************************

Aa= 0
Ba= 0
h1a= 0
h2a= 0

!*****************************************************************
!**********************Minimising Lamda1 (h1)*********************
deltah= 0.0001
h1b= 0.001
h1a= 0
DO 10
chisqa= 0
DO 20 i= 1,17
Cs(i)= exp(-h1a*t(i))!*Aa !+ Ba*exp(-h2a*t(i))
Ctdiff= Ct - Cs
c= Ctdiff**2/eCt**2
chisqa= chisqa + c(i)
20 END DO
!Use initial h1 value of 0 to calculate start-point chisq.

chisqb= 0
DO 30 i= 1,17
h1b= h1b + deltah
Cs(i)= exp(-h1b*t(i))!*Ab !+ Bb*exp(-h2b*t(i))
Ctdiff(i)= Ct(i) - Cs(i)
c= Ctdiff(i)**2/eCt(i)**2
chisqb= chisqb + c(i)
30 END DO
!First-step h1 used to find competing chisq for comparison.
WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1a
WRITE(6,*) 'Chi-squared b=', chisqb,'for Lamda1=',h1b
!Prints the two calculated chisq values to screen.
dchisq= chisqa - chisqb
IF (dchisq.GT.0) THEN
h1a= h1b
ELSE IF (dchisq.LE.0) THEN
deltah= deltah-((deltah*2)/100)

END IF
IF (chisqb.LE.6618.681) EXIT
10 END DO
WRITE(6,*) 'Chi-squared is', chisqb,'for Lamda1=',h1b
!*****************************************************************

END PROGRAM assignmentThanks for the help.
Nathan

EDIT: I should be getting a chi-squared of 6618.681 from this, but it's just stuck between 6921.866 and 6920.031. Help!

 

[Mark44]

Hey guys! Having issues with a code I'm writing. It doesn't fit the typical format and I may have already accidentally posted it in the wrong forum, so here's a copy of what was there with the newest version of the code. It'll (eventually) be a four-parameter least squares minimisation, but for the moment I'm just trying to get the single parameter version to work. I've got the first module done but, when I run it it gets caught up in a loop and doesn't actually minimise anything! Good start huh?

Here's my code, bless me with your insight!

      PROGRAM assignment
      
      ! A program designed to fit experiemental data, using the method
      ! of least squares to minimise the associated chi-squared and
      ! obtain the four control parameters A,B,h1 and h2.
      !*****************************************************************

      IMPLICIT NONE
      INTEGER i,t0
      DOUBLE PRECISION t(17),t2(17),Ct(17),eCt(17),Ctdiff(17),c(17)
      DOUBLE PRECISION Aa,Ab,Ba,Bb,h1a,h1b,h2a,h2b,chisqa,chisqb,dchisq
      DOUBLE PRECISION deltah,Cs(17)

        OPEN(21, FILE='data.txt', FORM='FORMATTED', STATUS='OLD')
           DO i=1,17
              READ(21,*)t(i),Ct(i),eCt(i)
           END DO
        CLOSE(21)

      !Read in data.txt as three one dimensional arrays.
      !*****************************************************************
      !OPEN(21, FILE='outtest.txt', FORM='FORMATTED', STATUS='NEW')
      !DO i=1,17
      !   WRITE(21,*)t(i),Ct(i),eCt(i)
      !END DO
      !CLOSE(21)
      !
      !Just to check input file is being read correctly.
      !*****************************************************************
      !****************************Parameters***************************
      
      Aa= 0
      Ba= 0
      h1a= 0
      h2a= 0

      !*****************************************************************
      !**********************Minimising Lamda1 (h1)*********************
      deltah= 0.0001
      h1b= 0.001
      h1a= 0
           DO 10
      chisqa= 0
              DO 20 i= 1,17   !  *** Is part of the line below supposed to be commented out?
                   Cs(i)= exp(-h1a*t(i))!*Aa !+ Ba*exp(-h2a*t(i))
                     Ctdiff= Ct - Cs  ! *** This is wrong. All three variables are arrays.
                       c= Ctdiff**2/eCt**2
                         chisqa= chisqa + c(i)
 20           END DO
      !Use initial h1 value of 0 to calculate start-point chisq.

      chisqb= 0
              DO 30 i= 1,17
                   h1b= h1b + deltah   ! *** In the line below, are you commenting out part of the line?
                     Cs(i)= exp(-h1b*t(i))!*Ab !+ Bb*exp(-h2b*t(i))
                       Ctdiff(i)= Ct(i) - Cs(i)
                         c= Ctdiff(i)**2/eCt(i)**2
                           chisqb= chisqb + c(i)
 30           END DO
      !First-step h1 used to find competing chisq for comparison.
                 WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1a
                 WRITE(6,*) 'Chi-squared b=', chisqb,'for Lamda1=',h1b
      !Prints the two calculated chisq values to screen.
                 dchisq= chisqa - chisqb
              IF (dchisq.GT.0) THEN
                 h1a= h1b
                 ELSE IF (dchisq.LE.0) THEN
                 deltah= deltah-((deltah*2)/100)

              END IF
           IF (chisqb.LE.6618.681) EXIT
 10        END DO
           WRITE(6,*) 'Chi-squared is', chisqb,'for Lamda1=',h1b
      !*****************************************************************
           
      END PROGRAM assignment

Thanks for the help.
Nathan

EDIT: I should be getting a chi-squared of 6618.681 from this, but it's just stuck between 6921.866 and 6920.031. Help!

I added some comments where I think you're going wrong. My comments look like this: ! ***
In one comment, you have three assignment statements where you are using array variables, but without indexes. This is wrong. Everywhere you have an array variable in a calculation, it should have an index. E.g., it should be Ct(i), not just Ct.

In my other comments, you have a line of code that might be wrong. It looks like you commented it out, but I'm not sure what you did.

Your variable names are terrible! For example, Aa,Ab,Ba,Bb,h1a,h1b. Can't you come up with some names that are self-explanatory? Also, why does your code meander all across the page? This makes it harder to see which statements are inside IF blocks and DO loops and such.