# Fortran77 woes- single parameter least squares minimisation stuck!

1. Jan 26, 2014

### NathanM

Hey guys! Having issues with a code I'm writing. 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!

Code (Text):
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
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
DO 10
chisqa= 0
DO 20 i= 1,17
h1a= h1a + deltah
Cs(i)= exp(-h1a*t(i))!*Aa !+ Ba*exp(-h2a*t(i))
Ctdiff= Ct - Cs
c= Ctdiff**2/eCt**2
chisqa= chisqa + c(i)
deltah= 0.001
20           END DO
!Use initial h1 value of 0 to calculate start-point chisq.

h1b= 0.001
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*2)/100

END IF
IF (chisqb.LE.6000) 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!

EDIT: I've read this is the wrong place for this. Sorry, have opened a thread in the correct place, please remove this one.

Last edited: Jan 26, 2014
2. Jan 26, 2014

### gsal

Code (Text):

1  PROGRAM assignment
2
3  ! A program designed to fit experiemental data, using the method
4  ! of least squares to minimise the associated chi-squared and
5  ! obtain the four control parameters A,B,h1 and h2.
6  !*****************************************************************
7
8  IMPLICIT NONE
9  INTEGER i,t0
10  DOUBLE PRECISION t(17),t2(17),Ct(17),eCt(17),Ctdiff(17),c(17)
11  DOUBLE PRECISION Aa,Ab,Ba,Bb,h1a,h1b,h2a,h2b,chisqa,chisqb,dchisq
12  DOUBLE PRECISION deltah,Cs(17)
13
14  OPEN(21, FILE='data.txt', FORM='FORMATTED', STATUS='OLD')
15  DO i=1,17
17  END DO
18  CLOSE(21)
19
20  !Read in data.txt as three one dimensional arrays.
21  !*****************************************************************
22  !OPEN(21, FILE='outtest.txt', FORM='FORMATTED', STATUS='NEW')
23  !DO i=1,17
24  ! WRITE(21,*)t(i),Ct(i),eCt(i)
25  !END DO
26  !CLOSE(21)
27  !
28  !Just to check input file is being read correctly.
29  !*****************************************************************
30  !****************************Parameters***************************
31
32  Aa= 0
33  Ba= 0
34  h1a= 0
35  h2a= 0
36
37  !*****************************************************************
38  !**********************Minimising Lamda1 (h1)*********************
39  deltah= 0
40  DO 10
41      chisqa= 0
42      DO 20 i= 1,17
43          h1a= h1a + deltah
44          Cs(i)= exp(-h1a*t(i))!*Aa !+ Ba*exp(-h2a*t(i))
45          Ctdiff= Ct - Cs
46          c= Ctdiff**2/eCt**2
47          chisqa= chisqa + c(i)
48          deltah= 0.001
49      20 END DO
50
51      !Use initial h1 value of 0 to calculate start-point chisq.
52      h1b= 0.001
53      chisqb= 0
54      DO 30 i= 1,17
55          h1b= h1b + deltah
56          Cs(i)= exp(-h1b*t(i))!*Ab !+ Bb*exp(-h2b*t(i))
57          Ctdiff(i)= Ct(i) - Cs(i)
58          c= Ctdiff(i)**2/eCt(i)**2
59          chisqb= chisqb + c(i)
60      30 END DO
61
62      !First-step h1 used to find competing chisq for comparison.
63      WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1a
64      WRITE(6,*) 'Chi-squared b=', chisqb,'for Lamda1=',h1b
65      !Prints the two calculated chisq values to screen.
66      dchisq= chisqa - chisqb
67      IF (dchisq.GT.0) THEN
68          h1a= h1b
69      ELSE IF (dchisq.LE.0) THEN
70          deltah= (-deltah*2)/100
71      END IF
72      IF (chisqb.LE.6000) EXIT
73  10 END DO
74  WRITE(6,*) 'Chi-squared is', chisqb,'for Lamda1=',h1b
75  !*****************************************************************
76
77  END PROGRAM assignment
78

First question...are you compiling with Fortran 77 or Fortran 90? I think 90.

I say this because you have typos in the program that may be doing something you did not intend. Recall that Fortran 90 operates on entire arrays when instructed to do so, i.e., no index is passed.

Using numbered linsting above:

line 44 calculates i-th index of Cs and only this index
line 45 performs an array operation!...is this what you want?
line 46 performs an array operation.
line 47 takes a single index from c

what is the point of performing arrays operations if ultimately you seem to just need a single index?

deltah is set to 0.001 at the end of do-loop 20

line 52 sets h1b to 0.001 and
line 53 sets chisqb to 0

this means that do-loop 30 does exactly the same thing on every iteration!

in do-loop 30, you seem to correctly use parenthesis and the i index in your operations, yet, in
line 58 you assign a single scalar to the entire 'c' array (since you did not specify an index)
line 59 will pick up the intended value and assign to chisqb

line 70 may flip the sign of deltah and divide it by 50, but this new value will be used only once at the beginning of do-loop 20; is this what you want? a short step back? just asking

3. Jan 26, 2014

### NathanM

Hi! Thanks for the reply :) I'm compiling with a program called Force2.0, but it's definitely meant to be fortran 77. With skills like mine, wouldn't be surprised if I've screwed up a few times!

I'll bullet point responses to keep things clear(ish).
(44)- Riiiight. So, having i=1,17 there isn't running through all three equations and then summing? Crap.
(45,46,47)- Yeah I think so- this whole do block is supposed to be calculating chi-squared for values of t(i)=1,17 (ie all seventeen values of t read in from the data file), then adding them all together to find the actual chi-squared for that particular value of lamda 1 (h1 as I've entered it).

NB: I can't believe I didn't notice that the middle two equations are missing their array assignments. Dumbass. All four should be using arrays.

- Yeah I wasn't sure of the placement of deltah at the time of posting, I've moved it outside of the main loop since then. I'll post my newest code at the end of this reply. Likewise with h1b- that's been moved out the way too. I've left chisqb there as it makes sense in my head that chisq should be zeroed each time.

- Another delightful cockup on my part, thanks for pointing it out. These missing arrays must be from when I imported them from my original attempt (tried using functions and subroutines, did not end well), but still can't believe I (nor my friends who had a quick look) didn't notice them.

- hmm no. The idea was it's supposed to run through the whole loop (only operates in the second inner loop in the new code), knock the step back and downsize, but only for that iteration (unless dchisq is negative again). I've edited it to take two steps back now and decrease step size further (though I may change the size of the decrease once the minimisation loop is working).

Here's the newest code (I've added back in the missing array identifiers).

Code (Text):
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
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
h1a= h1a
Cs(i)= exp(-h1a*t(i))!*Aa !+ Ba*exp(-h2a*t(i))
Ctdiff(i)= Ct(i) - Cs(i)
c(i)= Ctdiff(i)**2/eCt(i)**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(i)= 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
Thankyou for the help :)

EDIT: I imagine I'm going to have to use a function for the chi-squared calculation aren't I?

Last edited: Jan 26, 2014
4. Jan 26, 2014

### Staff: Mentor

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.

5. Jan 26, 2014

### NathanM

Hi! I'll reply using bullet points for ease of reading. As always, thanks for taking the time to read through my mess :)

- RE: commented out variables- yeah I've commented them out whilst I work on the single parameter minimisation. They'll come back in one by one once it's working.
- You caught my cockup too. No clue how I managed to go the last few hours without noticing. Guess I just didn't expect there to be issues with those blocks and ignored them- foolish.

As for variable names, well, they made sense to me at the time but if you'd prefer something clearer I can alter them. RE: meandering code, I assume you mean the equation blocks? They weren't like that before, got whinged at by some here when they had a look so changed it, but I'll change it back.
If you mean the code in general, I thought it was okay. Will alter it if you have any specific preferences.

Last edited: Jan 26, 2014
6. Jan 26, 2014

### Staff: Mentor

When you post code, it's good practice here at PF to put [code] ... [/code] tags around it. I've done this for you in your code.
Probably be a good idea.
Is your code producing the right results now?
Don't do it for me - do it for yourself. When you review this code in a month or two, having meaningless variable names will make it harder for you to understand the code, even though you wrote it yourself.
We call them assignments statements.
Minor nit, but this Greek letter - λ - is written "lambda".

7. Jan 26, 2014

### NathanM

Done!

So as the assignment statements stand right now, are they not running through the calculation for every array labelled parameter using the i=1,17 statement? If no, it'll be better to move all of the assignment statements to a function and just execute that 17 times, yes? I'm not great with functions; going to throw one together quick and post it in an edit here, help with any needed corrections would be appreciated.

No, still stuck at two lamda values.

This is true. What might be a great idea at the time might not be later. Not particularly creative, but I've just renamed them after the do loop they're used in.
I know re: lamda, but I didn't want "lamda1, lamda2" lengthening lines out. I usually use h in its place in fortran, since the actual character isn't available.

8. Jan 26, 2014

### NathanM

I've written this function to replace the assignment statements in the two interior loops, leaving the summation statements in.

Code (Text):

DO 20 i= 1,17
h1a= h1a
chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),Cs(i),eCt(i))
20           END DO

DOUBLE PRECISION FUNCTION c(Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt)
IMPLICIT NONE
DOUBLE PRECISION Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt

Cs= exp(-h1loop1*t)!*Aa !+ Ba*exp(-h2loop1*t(i))
Ctdiff= Ct - Cs
c= Ctdiff**2/eCt**2

RETURN
END
Will this do the trick?

EDIT: It doesn't like the array indexes in the function call-in in the main program body. Not quite sure how to properly call it in...do I just assign the indexes inside the function itself (swear somebody told me you couldn't use them in functions) and write the call-in as c(Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt)?

Last edited: Jan 26, 2014
9. Jan 27, 2014

### Staff: Mentor

There are several things wrong with what you're doing. The DO loop is presumably part of your main program. It would be in some code that looks like this:
Code (Text):

PROGRAM Whatever

! Declarations
...
DO 20 i= 1,17
h1a= h1a
chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),Cs(i),eCt(i))
20 END DO
...
END PROGRAM Whatever

In the code above you have a completely useless assignment statement: h1a = h1a. Here you are assigning the value of h1a to itself. That value is already there, so there's no point in doing this.

In your main program you need a declaration for your c function. BTW, you really should have a better name for it than 'c'. For a function, a reasonable name is a noun that describes the value that is being returned.

Below the main program are any subprograms, including your function, which is one type of subprogram (subroutine being the other).

The main program has the statement that is the call to the function. (In 40+ years I've never heard anybody refer to it as a 'call-in'.) Your function call looks OK - it has the right number of parameters, and they appear to be the correct types. After you have your function declaration in place, I think things should work properly.

It's lamBda - there's a 'B' in the name.

10. Jan 27, 2014

### NathanM

I wasn't sure if h1a would need stating before the loop or not, now I do.
c works for me. It's an intermediate value that has no use outside of its function, and it's only definition is as the sum of c(i) array- Chi-squared.
I like to be different. That's what it does so works for me to call it that, I'm a physicist first not a programmer.

Here's the code as it now stands.

Code (Text):
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)
DOUBLE PRECISION Aloop1,Aloop2,Bloop1,Bloop2,h1loop1,h1loop2
DOUBLE PRECISION deltah,Cs(17),h2loop1,h2loop2
DOUBLE PRECISION chisqa,chisqb,dchisq,c(17,17,17,17,17,17,17)

OPEN(21, FILE='data.txt', FORM='FORMATTED', STATUS='OLD')
DO i=1,17
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***************************

Aloop1= 0
Bloop2= 0
h1loop1= 0
h2loop2= 0

!*****************************************************************
!**********************Minimising Lamda1 (h1)*********************
deltah= 0.0001
h1loop2= 0.001
h1loop1= 0
DO 10
chisqa= 0
DO 20 i= 1,17
chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
>           Cs(i),eCt(i))
20           END DO
WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1loop1
!Use initial h1 value of 0 to calculate start-point chisq.

chisqb= 0
DO 30 i= 1,17
h1loop2= h1loop2 + deltah
chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
>           Cs(i),eCt(i))
30           END DO
!First-step h1 used to find competing chisq for comparison.
WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1loop1
WRITE(6,*) 'Chi-squared b=', chisqb,'for Lamda1=',h1loop2
!Prints the two calculated chisq values to screen.
dchisq= chisqa - chisqb
IF (dchisq.GT.0) THEN
h1loop1= h1loop2
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=',h1loop2
!*****************************************************************

END PROGRAM assignment

!*****************************************************************
!***************************FUNCTIONS*****************************
DOUBLE PRECISION FUNCTION c(Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt)
IMPLICIT NONE
DOUBLE PRECISION Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt

Cs= exp(-h1loop1*t)!*Aa !+ Ba*exp(-h2loop1*t(i))
Ctdiff= Ct - Cs
c= Ctdiff**2/eCt**2

RETURN
END FUNCTION c(Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt)
Just took a stab at defining c at the start of the program, wasn't entirely sure but there you go. This presents the following errors on compiling.

test.f:46.36:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:46.42:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:46.50:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
C:\Users\Nathan\Desktop\Documents\University of Birmingham MSc Physics and Technology of Nuclear Reactors\Lectures\fortran\test.f:46.55:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:46.65:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:46.71:

chisqa= chisqa + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:47.23:

> Cs(i),eCt(i))
1
Warning: Extension: REAL array index at (1)
test.f:55.36:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:55.42:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:55.50:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:55.55:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:55.65:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:55.71:

chisqb= chisqb + c(Cs(i),h1loop1,t(i),Ctdiff(i),Ct(i),
1
Warning: Extension: REAL array index at (1)
test.f:56.23:

> Cs(i),eCt(i))
1
Warning: Extension: REAL array index at (1)
test.f:1.72:

PROGRAM assignment
1
test.f:78.72:

DOUBLE PRECISION FUNCTION c(Cs,h1loop1,t,Ctdiff,Ct,Cs,eCt)
2
Error: Two main PROGRAMs at (1) and (2)
test.f:13: error: size of variable 'c' is too large

Surprised by the last one. I thought fortran 77 could handle seven dimensions to an array. Though perhaps I've just not entered something correctly.

11. Jan 27, 2014

### gsal

O.k., stay away from functions, for now. But please, just because you are getting help here, don't stop studying Fortran...you really need to for this assignment and as a physicist!

It seems that your comments in the source apply to the code above them...that is not how it is typically done; instead, you first tell the reader what you are about to do, then you do it.

From what I can tell, the only variables that need to be arrays are the ones you read from the file; so, I cleaned that.

Deleted unnecessary stuff to force you to add it later with consistent variable naming.

Changed the program a bit to force you to re-read it and focus on the correctness, once again.

And that is it, I am not going to bother to explain EVERYTHING that is wrong with your latest code (you really went backwards); instead, please take a good look at the cleaned-up code that I post below...it at least compiles; but, as simple as it might be, I have not check the math at all or the looping or anything, I just looked at the Fortran.

Code (Text):

!
! 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.
!
PROGRAM assignment

IMPLICIT NONE
INTEGER i
DOUBLE PRECISION t(17), Ct(17), eCt(17)
DOUBLE PRECISION h1loop1, h1loop2, deltah, Cs
DOUBLE PRECISION chisqa, chisqb, dchisq

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

!**********************Minimising Lamda1 (h1)*********************
deltah  = 0.0001
h1loop2 = 0.001
h1loop1 = 0.0 !Use initial value of 0 to calculate start-point chisq
DO 10
chisqa = 0.0
DO 20 i = 1, 17
Cs     = exp(-h1loop1*t(i))
chisqa = chisqa + ( (Ct(i) - Cs)/eCt(i) )**2
20       END DO

chisqb = 0.0
DO 30 i = 1, 17
h1loop2 = h1loop2 + deltah
Cs      = exp(-h1loop2*t(i))
chisqb  = chisqb + ( (Ct(i) - Cs)/eCt(i) )**2
30       END DO

!Print the two calculated chisq values to screen.
WRITE(6,*) 'Chi-squared a=', chisqa,'for Lamda1=',h1loop1
WRITE(6,*) 'Chi-squared b=', chisqb,'for Lamda1=',h1loop2

dchisq = chisqa - chisqb
IF (dchisq.GT.0.0) THEN
h1loop1 = h1loop2
ELSE
deltah = deltah-((deltah*2)/100)
END IF

IF (chisqb.LE.6618.681) EXIT
10   END DO

WRITE(6,*) 'Chi-squared is', chisqb,' for Lamda1 = ', h1loop2
END PROGRAM assignment

12. Jan 27, 2014

### NathanM

lol oh don't worry, I'm going to be getting a lot more familiar with fortran77 over the next few years. We're only really getting into it now, but it'll be used in my lab reports later in the semester and my dissertation over the summer. Then, I hear reactor control is still done in fortran77...there'll be no escape!

No probs, I'll move the comments around. I'll compare the new code to my own latest and learn from it, thanks. Updated my code, the iteration is still getting stuck between two values. This problem is ridiculous.

Last edited: Jan 27, 2014
13. Jan 27, 2014

### SteamKing

Staff Emeritus
You have defined 'c' as a DOUBLE PRECISION array in the main program and then again as a DOUBLE PRECISION function. Which is it? That's why you are getting a bunch of compiler errors.

14. Jan 27, 2014

### gsal

I would recommend reducing the number of points from 17 down to 3 and performing a hand calculation...I mean, really, hand calculate values one iteration at a time to really see what your equations are yielding, in particular, when to increment h1loop1 and h1loop2.