1c program DRSNLSGB 2c>> 2001-05-24 DRSNLSGB Krogh Minor change for making .f90 version. 3c>> 1997-06-18 DRSNLSGB Krogh Changes to improve C portability. 4c>> 1994-11-02 DRSNLSGB Krogh Changes to use M77CON 5c>> 1994-09-14 DRSNLSGB CLL Set IV(OUTLEV) = 0 for comparing output. 6c>> 1992-04-13 CLL Rename and reorder common block [D/S]KEY. 7c>> 1992-02-03 CLL @ JPL 8c>> 1990-07-02 CLL @ JPL 9c>> 1990-06-27 CLL @ JPL 10c>> 1990-04-05 CLL @ JPL 11c>> 1990-03-29 CLL @ JPL 12c Demo driver for SNLSGB. A variant of the nonlinear LS code NL2SOL. 13c SNLSGB solves the "separable" problem. 14c SNLSGB handles bounds on the nonlinear variables. 15c SNLSGB requires values of the function and the Jacobian matrix. 16c Note: MDER is the number of ones in the array IND(). 17c ------------------------------------------------------------------ 18c--S replaces "?": DR?NLSGB, ?NLSGB, ?CALCA, ?CALCB, ?IVSET, ?KEY 19c ------------------------------------------------------------------ 20 external SCALCA, SCALCB 21 integer J,IVAR,ITERM,LIV,LV,MA,MB,MDATA,MDER,NA,NB,NDATA 22 parameter(MDATA = 30, MA = 2, MB = 5) 23 parameter(MDER = 4) 24 parameter(LIV = 115 + 4*MA + MB + 2*MDER) 25 parameter(LV = 105 + MDATA*(MB+MDER+MA+3) + 26 * (MB*(MB+3))/2 + MA*(2*MA+21)) 27 integer IND(MB+1,MA), IV(LIV) 28 integer F, OUTLEV, SOLPRT, STATPR, X0PRT 29 parameter(F=10, OUTLEV=19, SOLPRT=22, STATPR=23, X0PRT=24) 30 real ALF(MA), BET(MB), BND(2,MA) 31 real DOF, V(LV), YDATA(MDATA) 32 data (BND(1,J),J=1,2) / 2*5.0e0 / 33 data (BND(2,J),J=1,2) / 2*10.0e0 / 34 data ((IND(ITERM,IVAR),IVAR = 1,2),ITERM = 1,6)/ 35 * 0, 0, 36 * 1, 0, 37 * 1, 0, 38 * 0, 1, 39 * 0, 1, 40 * 0, 0/ 41 data YDATA / 42 * 1.700641e0, 1.793512e0, 1.838309e0, 1.838416e0, 1.792204e0, 43 * 1.700501e0, 1.579804e0, 1.426268e0, 1.260724e0, 1.084901e0, 44 * 0.917094e0, 0.761920e0, 0.627304e0, 0.522146e0, 0.446645e0, 45 * 0.404920e0, 0.392033e0, 0.409622e0, 0.453045e0, 0.510765e0, 46 * 0.584554e0, 0.663109e0, 0.747613e0, 0.829439e0, 0.908496e0, 47 * 0.983178e0, 1.051046e0, 1.114072e0, 1.171746e0, 1.227823e0/ 48c ------------------------------------------------------------------ 49 NDATA = MDATA 50 NA = MA 51 NB = MB 52 ALF(1) = 5.0e0 53 ALF(2) = 10.0e0 54 IV(1) = 0 55 56 print'(1x,a)', 57 * 'Program DRSNLSGB.. Demo driver for SNLSGB.', 58 * ' A variant of NL2SOL.', 59 * ' SNLSGB handles the Separable problem.', 60 * ' SNLSGB requires values of the function and the Jacobian.', 61 * ' SNLSGB handles bounds on the nonlinear variables.', 62 * ' ', 63 * 'Sample problem is a nonlinear curve fit to data.', 64 * 'Model function is B1 + B2 * cos(A1*t) + B3 * sin(A1*t) +', 65 * ' B4 * cos(A2*t) + B5 * sin(A2*t) + Noise', 66 * 'Data generated using', 67 * '(A1, A2, B1, ..., B5) = (6, 9, 1, 0.5, 0.4, 0.2, 0.1)', 68 * 'and Gaussian noise with mean 0 and', 69 * 'sample standard deviation 0.001', 70 * ' ' 71 72 call SIVSET(1, IV, LIV, LV, V) 73 IV( X0PRT) = 1 74 IV(OUTLEV) = 0 75 IV(STATPR) = 1 76 IV(SOLPRT) = 1 77 78 call SNLSGB(NDATA, NA, NB, ALF, BND, BET, YDATA, SCALCA, SCALCB, 79 * IND, NB+1, IV, LIV, LV, V) 80 81 DOF = max(NDATA - NA - NB, 1) 82 print'(1x/1x,a,g12.4)', 83 * 'SIGFAC: sqrt((2 * V(F))/DOF) =', 84 * sqrt(2.0e0 * V(F)/DOF) 85 stop 86 end 87c ================================================================== 88 subroutine SCALCA(NDATA, NA, NB, ALF, NCOUNT, PHI) 89c Test case for separable nonlinear least squares computation. 90c Computes MDATA x NB matrix PHI as a function of the 91c nonlinear parameters ALF(). 92c For J .le. NB the (I,J) term of PHI is the coefficient of the 93c linear coefficient B(J) in row I of the model. 94c In this example the model does not have a term that is not 95c multiplied by a linear coefficient. If such a term is present 96c then PHI must have an (NB+1)st column to hold this term. 97c This code Illustrates saving results in common between SCALCA and 98c SCALCB to avoid recalculation of common subexpressions. 99c ------------------------------------------------------------------ 100 common/SKEY/C1, C2, S1, S2, KEY 101 save /SKEY/ 102 integer I, KEY, MDATA, NA, NB, NCOUNT, NDATA 103 parameter(MDATA = 30) 104 real ALF(NA), C1(MDATA), C2(MDATA) 105 real DEL, PHI(NDATA,NB) 106 real S1(MDATA), S2(MDATA), T 107c ------------------------------------------------------------------ 108 T = 0.0E0 109 DEL = 1.0E0 / 29.0E0 110 KEY = NCOUNT 111 do 10 I = 1,NDATA 112 C1(I) = cos(ALF(1)*T) 113 S1(I) = sin(ALF(1)*T) 114 C2(I) = cos(ALF(2)*T) 115 S2(I) = sin(ALF(2)*T) 116 PHI(I,1) = 1.0E0 117 PHI(I,2) = C1(I) 118 PHI(I,3) = S1(I) 119 PHI(I,4) = C2(I) 120 PHI(I,5) = S2(I) 121 T = T + DEL 122 10 continue 123 return 124 end 125c ================================================================== 126 subroutine SCALCB(NDATA, NA, NB, ALF, NCOUNT, DER) 127c Test case for separable nonlinear least squares computation. 128c Computes the NDATA x NDER matrix DER. Here NDER is the number of 129c ones in the indicator array IND(). The columns of DER correspond 130c to nonzero entries of IND() traversed columnwise. 131c In this example NDER is 4 and the correspondence is as follows: 132c Col of DER: 1 2 3 4 133c Element of IND(): (2,1) (3,1) (4,2) (5,2) 134c In this example Row I of DER will be set to contain the values 135c of the four partial derivatives: 136c Partial of PHI(I,2) with respect to ALP(1) 137c Partial of PHI(I,3) with respect to ALP(1) 138c Partial of PHI(I,4) with respect to ALP(2) 139c Partial of PHI(I,5) with respect to ALP(2) 140c ------------------------------------------------------------------ 141 common/SKEY/C1, C2, S1, S2, KEY 142 save /SKEY/ 143 integer I, KEY, MDATA, NCOUNT, NDATA, NA, NB, NDER 144 parameter(MDATA = 30, NDER = 4) 145 real ALF(NA), C1(MDATA), C2(MDATA) 146 real DEL, DER(NDATA,NDER), S1(MDATA), S2(MDATA), T 147c ------------------------------------------------------------------ 148 T = 0.0E0 149 DEL = 1.0E0 / 29.0E0 150 if(NCOUNT .eq. KEY) then 151 do 10 I = 1,NDATA 152 DER(I,1) = -S1(I)*T 153 DER(I,2) = C1(I)*T 154 DER(I,3) = -S2(I)*T 155 DER(I,4) = C2(I)*T 156 T = T + DEL 157 10 continue 158 else 159 do 20 I = 1,NDATA 160 DER(I,1) = -sin(ALF(1)*T)*T 161 DER(I,2) = cos(ALF(1)*T)*T 162 DER(I,3) = -sin(ALF(2)*T)*T 163 DER(I,4) = cos(ALF(2)*T)*T 164 T = T + DEL 165 20 continue 166 endif 167 return 168 end 169