V3FIT
lmpar_parallel.f90
1  SUBROUTINE lmpar_parallel(j, fcn)
2  USE fdjac_mod, ONLY: wa1p => wap, m=>mp, n=>np, &
3  num_lm_params, xp, ncnt=>ncntp
4  USE lmpar_mod
5  IMPLICIT NONE
6 !-----------------------------------------------
7 ! D u m m y A r g u m e n t s
8 !-----------------------------------------------
9  INTEGER :: j
10  EXTERNAL fcn
11 #if !defined(MPI_OPT)
12 !-----------------------------------------------
13 ! L o c a l P a r a m e t e r s
14 !-----------------------------------------------
15  REAL(rprec), PARAMETER :: zero = 0
16  REAL(rprec), DIMENSION(11), PARAMETER :: factors = &
17  (/ 1.0_dp, 0.5_dp, 0.25_dp, 0.128_dp, 0.75_dp, &
18  1.25_dp, 1.5_dp, 0.9_dp, 1.1_dp, 1.75_dp, 2.1_dp /)
19 !-----------------------------------------------
20 ! L o c a l V a r i a b l e s
21 !-----------------------------------------------
22  INTEGER :: iflag, nfact
23 #if defined(CRAY)
24  INTEGER :: istat, k
25 #endif
26  REAL(rprec) :: deltain, parin, fnorm_in, pnorm_in, scale_factor
27 !-----------------------------------------------
28 ! E x t e r n a l F u n c t i o n s
29 !-----------------------------------------------
30  REAL(rprec), EXTERNAL :: enorm
31 !-----------------------------------------------
32 !
33 ! THIS ROUTINE IS PASSED TO THE MULTI-PROCESSOR HANDLING
34 ! ROUTINE
35 
36 
37 ! ***************************************************
38 ! stepping algorithm similar to that used in the original parallel optimizer
39 ! by M.Zarnstorff and S. Ethier, Feb. 1999
40 !
41 ! Re-implemented, MCZ July 2000
42 ! ***************************************************
43  nfact = SIZE(factors)
44 
45  IF (lfirst_lm .and. num_lm_params > 2) THEN
46 !
47 ! do an exponential spread the first time to see where we are
48 !
49 !SPH scale_factor = EXP((j-1)*log(spread_ratio)/num_lm_params)
50  scale_factor = 10._dp**(1._dp - j)
51  ELSE IF (num_lm_params > 2*nfact) THEN
52  scale_factor = (j*maxval(factors))/num_lm_params
53  ELSE IF (j .le. nfact) THEN
54  scale_factor = factors(j)
55  ELSE
56  scale_factor =((j-nfact)*minval(factors))/(num_lm_params-nfact)
57  ENDif
58 
59  deltain = delta * scale_factor
60 
61 !
62 ! Compute perturbation vector (wa1p) and Lev/Marq PARAMETER (par)
63 ! for different tolerances, delta
64 !
65  parin = par
66 
67  CALL lmpar (n, fjac, ldfjac, ipvt, diag, qtf, deltain, parin, &
68  wa1p, wa2p, wa3p, wa4p)
69 
70 !
71 ! store the direction p and x + p. calculate the norm of p.
72 !
73  IF (parin.eq.zero .and. j.ne.1) wa1p = wa1p*scale_factor
74  wa1p = -wa1p
75  wa2p = xp + wa1p
76  wa3p = diag*wa1p
77  pnorm_in = enorm(n, wa3p)
78 
79 !
80 ! evaluate the function at x + p and calculate its norm.
81 !c
82  iflag = j
83  CALL fcn (m, n, wa2p, wa4p, iflag, ncnt)
84 
85  fnorm_in = enorm(m, wa4p)
86 
87 !
88 ! OPEN A UNIQUE FILE FOR I/O IN MULTI-PROCESSOR SYSTEM
89 !
90  WRITE (j+1000) j, iflag, pnorm_in, fnorm_in, parin, deltain
91 #if defined(CRAY)
92  DO k = 1, n
93  WRITE (j+1000) wa1p(k), wa2p(k)
94  DO istat = 1, n
95  WRITE (j+1000) fjac(k, istat)
96  END DO
97  END DO
98  DO k = 1, m
99  WRITE (j+1000) wa4p(k)
100  END DO
101 #else
102  WRITE (j+1000) wa1p, wa2p, wa4p, fjac(1:n, 1:n)
103 #endif
104  CLOSE (j+1000) !!Needed to run correctly in multi-tasking...
105 #endif
106  END SUBROUTINE lmpar_parallel