source: LMDZ6/branches/Amaury_dev/libf/phylmd/read_pstoke0.F90 @ 5103

Last change on this file since 5103 was 5103, checked in by abarral, 8 weeks ago

Handle CPP_INLANDSIS in lmdz_cppkeys_wrapper.F90
Remove obsolete key wrgrads_thermcell, _ADV_HALO, _ADV_HALLO, isminmax
Remove redundant uses of CPPKEY_INCA (thanks acozic)
Remove obsolete misc/write_field.F90
Remove unused ioipsl_* wrappers
Remove calls to WriteField_u with wrong signature
Convert .F -> .[fF]90
(lint) uppercase fortran operators
[note: 1d and iso still broken - working on it]
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 15.0 KB
Line 
1
2! $Id: read_pstoke0.F90 5103 2024-07-23 13:29:36Z abarral $
3
4
5
6SUBROUTINE read_pstoke0(irec, zrec, zkon, zkev, airefi, phisfi, t, mfu, mfd, &
7    en_u, de_u, en_d, de_d, coefh, fm_therm, en_therm, frac_impa, frac_nucl, &
8    pyu1, pyv1, ftsol, psrf)
9
10  ! ******************************************************************************
11  ! Frederic HOURDIN, Abderrahmane IDELKADI
12  ! Lecture des parametres physique stockes online necessaires pour
13  ! recalculer offline le transport des traceurs sur la meme grille que
14  ! online
15  ! A FAIRE : une seule routine au lieu de 2 (lectflux, redecoupe)!
16  ! ******************************************************************************
17
18  USE netcdf, ONLY: nf90_inq_dimid,nf90_inquire_dimension,nf90_get_var,nf90_inq_varid,nf90_open,&
19          nf90_nowrite
20  USE dimphy
21  USE indice_sol_mod
22  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, nbp_lev
23
24  IMPLICIT NONE
25
26  INTEGER kon, kev, zkon, zkev
27!  PARAMETER (kon=iim*(jjm-1)+2, kev=llm)
28  REAL :: phisfi(nbp_lon*(nbp_lat-2)+2) !phisfi(kon)
29  REAL,ALLOCATABLE :: phisfi2(:,:) !phisfi2(nbp_lon, nbp_lat)
30  REAL,ALLOCATABLE :: airefi2(:,:) !airefi2(nbp_lon, nbp_lat)
31
32  REAL :: mfu(nbp_lon*(nbp_lat-2)+2,nbp_lev) !mfu(kon, kev)
33  REAL :: mfd(nbp_lon*(nbp_lat-2)+2,nbp_lev) !mfd(kon, kev)
34  REAL :: en_u(nbp_lon*(nbp_lat-2)+2,nbp_lev) !en_u(kon, kev)
35  REAL :: de_u(nbp_lon*(nbp_lat-2)+2,nbp_lev) !de_u(kon, kev)
36  REAL :: en_d(nbp_lon*(nbp_lat-2)+2,nbp_lev) !en_d(kon, kev)
37  REAL :: de_d(nbp_lon*(nbp_lat-2)+2,nbp_lev) !de_d(kon, kev)
38  REAL :: coefh(nbp_lon*(nbp_lat-2)+2,nbp_lev) !coefh(kon, kev)
39
40  ! abd 25 11 02
41  ! Thermiques
42  REAL :: fm_therm(nbp_lon*(nbp_lat-2)+2,nbp_lev) !fm_therm(kon, kev)
43  REAL :: en_therm(nbp_lon*(nbp_lat-2)+2,nbp_lev) !en_therm(kon, kev)
44  REAL :: t(nbp_lon*(nbp_lat-2)+2,nbp_lev) !t(kon, kev)
45
46  REAL,ALLOCATABLE :: mfu2(:,:,:) !mfu2(nbp_lon, nbp_lat, kev)
47  REAL,ALLOCATABLE :: mfd2(:,:,:) !mfd2(nbp_lon, nbp_lat, kev)
48  REAL,ALLOCATABLE :: en_u2(:,:,:) !en_u2(nbp_lon, nbp_lat, kev)
49  REAL,ALLOCATABLE :: de_u2(:,:,:) !de_u2(nbp_lon, nbp_lat, kev)
50  REAL,ALLOCATABLE :: en_d2(:,:,:) !en_d2(nbp_lon, nbp_lat, kev)
51  REAL,ALLOCATABLE :: de_d2(:,:,:) !de_d2(nbp_lon, nbp_lat, kev)
52  REAL,ALLOCATABLE :: coefh2(:,:,:) !coefh2(nbp_lon, nbp_lat, kev)
53  REAL,ALLOCATABLE :: t2(:,:,:) !t2(nbp_lon, nbp_lat, kev)
54  ! Thermiques
55  REAL,ALLOCATABLE :: fm_therm2(:,:,:) !fm_therm2(nbp_lon, nbp_lat, kev)
56  REAL,ALLOCATABLE :: en_therm2(:,:,:) !en_therm2(nbp_lon, nbp_lat, kev)
57
58  REAL,ALLOCATABLE :: pl(:) !pl(kev)
59  INTEGER irec
60  INTEGER xid, yid, zid, tid
61  INTEGER zrec, zim, zjm
62  INTEGER ncrec, nckon, nckev, ncim, ncjm
63
64  REAL :: airefi(nbp_lon*(nbp_lat-2)+2) !airefi(kon)
65  CHARACTER *20 namedim
66
67  ! !! attention !!
68  ! attention il y a aussi le pb de def kon
69  ! dim de phis??
70
71  REAL :: frac_impa(nbp_lon*(nbp_lat-2)+2,nbp_lev) !frac_impa(kon, kev)
72  REAL :: frac_nucl(nbp_lon*(nbp_lat-2)+2,nbp_lev) !frac_nucl(kon, kev)
73  REAL,ALLOCATABLE :: frac_impa2(:,:,:) !frac_impa2(nbp_lon, nbp_lat, kev)
74  REAL,ALLOCATABLE :: frac_nucl2(:,:,:) !frac_nucl2(nbp_lon, nbp_lat, kev)
75  REAL :: pyu1(nbp_lon*(nbp_lat-2)+2) !pyu1(kon)
76  REAL :: pyv1(nbp_lon*(nbp_lat-2)+2) !pyv1(kon)
77  REAL,ALLOCATABLE :: pyu12(:,:), pyv12(:,:) !pyu12(nbp_lon, nbp_lat), pyv12(nbp_lon, nbp_lat)
78  REAL :: ftsol(nbp_lon*(nbp_lat-2)+2,nbp_lev) !ftsol(kon, nbsrf)
79  REAL :: psrf(nbp_lon*(nbp_lat-2)+2,nbp_lev) !psrf(kon, nbsrf)
80  REAL,ALLOCATABLE :: ftsol1(:),ftsol2(:) !ftsol1(kon), ftsol2(kon)
81  REAL,ALLOCATABLE :: ftsol3(:),ftsol4(:) !ftsol3(kon), ftsol4(kon)
82  REAL,ALLOCATABLE :: psrf1(:), psrf2(:) !psrf1(kon), psrf2(kon)
83  REAL,ALLOCATABLE :: psrf3(:), psrf4(:) !psrf3(kon), psrf4(kon)
84  REAL,ALLOCATABLE :: ftsol12(:,:) !ftsol12(nbp_lon, nbp_lat)
85  REAL,ALLOCATABLE :: ftsol22(:,:) !ftsol22(nbp_lon, nbp_lat)
86  REAL,ALLOCATABLE :: ftsol32(:,:) !ftsol32(nbp_lon, nbp_lat)
87  REAL,ALLOCATABLE :: ftsol42(:,:) !ftsol42(nbp_lon, nbp_lat)
88  REAL,ALLOCATABLE :: psrf12(:,:) !psrf12(nbp_lon, nbp_lat)
89  REAL,ALLOCATABLE :: psrf22(:,:) !psrf22(nbp_lon, nbp_lat)
90  REAL,ALLOCATABLE :: psrf32(:,:) !psrf32(nbp_lon, nbp_lat)
91  REAL,ALLOCATABLE :: psrf42(:,:) !psrf42(nbp_lon, nbp_lat)
92
93  INTEGER,SAVE :: ncidp
94  INTEGER,SAVE :: varidmfu, varidmfd, varidps, varidenu, variddeu
95  INTEGER,SAVE :: varidt
96  INTEGER,SAVE :: varidend, varidded, varidch, varidfi, varidfn
97  ! therm
98  INTEGER,SAVE :: varidfmth, varidenth
99  INTEGER,SAVE :: varidyu1, varidyv1, varidpl, varidai, varididvt
100  INTEGER,SAVE :: varidfts1, varidfts2, varidfts3, varidfts4
101  INTEGER,SAVE :: varidpsr1, varidpsr2, varidpsr3, varidpsr4
102
103  INTEGER l, i
104  INTEGER start(4), count(4), status
105  REAL rcode
106  LOGICAL,SAVE :: first=.TRUE.
107
108  ! Allocate arrays
109  kon=nbp_lon*(nbp_lat-2)+2
110  kev=nbp_lev
111
112  ALLOCATE(phisfi2(nbp_lon, nbp_lat))
113  ALLOCATE(airefi2(nbp_lon, nbp_lat))
114  ALLOCATE(mfu2(nbp_lon, nbp_lat, kev))
115  ALLOCATE(mfd2(nbp_lon, nbp_lat, kev))
116  ALLOCATE(en_u2(nbp_lon, nbp_lat, kev))
117  ALLOCATE(de_u2(nbp_lon, nbp_lat, kev))
118  ALLOCATE(en_d2(nbp_lon, nbp_lat, kev))
119  ALLOCATE(de_d2(nbp_lon, nbp_lat, kev))
120  ALLOCATE(coefh2(nbp_lon, nbp_lat, kev))
121  ALLOCATE(t2(nbp_lon, nbp_lat, kev))
122  ALLOCATE(fm_therm2(nbp_lon, nbp_lat, kev))
123  ALLOCATE(en_therm2(nbp_lon, nbp_lat, kev))
124  ALLOCATE(pl(kev))
125  ALLOCATE(frac_impa2(nbp_lon, nbp_lat, kev))
126  ALLOCATE(frac_nucl2(nbp_lon, nbp_lat, kev))
127  ALLOCATE(pyu12(nbp_lon, nbp_lat), pyv12(nbp_lon, nbp_lat))
128  ALLOCATE(ftsol1(kon), ftsol2(kon))
129  ALLOCATE(ftsol3(kon), ftsol4(kon))
130  ALLOCATE(psrf1(kon), psrf2(kon))
131  ALLOCATE(psrf3(kon), psrf4(kon))
132  ALLOCATE(ftsol12(nbp_lon, nbp_lat))
133  ALLOCATE(ftsol22(nbp_lon, nbp_lat))
134  ALLOCATE(ftsol32(nbp_lon, nbp_lat))
135  ALLOCATE(ftsol42(nbp_lon, nbp_lat))
136  ALLOCATE(psrf12(nbp_lon, nbp_lat))
137  ALLOCATE(psrf22(nbp_lon, nbp_lat))
138  ALLOCATE(psrf32(nbp_lon, nbp_lat))
139  ALLOCATE(psrf42(nbp_lon, nbp_lat))
140
141  ! ---------------------------------------------
142  ! Initialisation de la lecture des fichiers
143  ! ---------------------------------------------
144
145  IF (irec==0) THEN
146
147    rcode = nf90_open('phystoke.nc', nf90_nowrite, ncidp)
148
149    rcode = nf90_inq_varid(ncidp, 'phis', varidps)
150    PRINT *, 'ncidp,varidps', ncidp, varidps
151
152    rcode = nf90_inq_varid(ncidp, 'sig_s', varidpl)
153    PRINT *, 'ncidp,varidpl', ncidp, varidpl
154
155    rcode = nf90_inq_varid(ncidp, 'aire', varidai)
156    PRINT *, 'ncidp,varidai', ncidp, varidai
157
158    rcode = nf90_inq_varid(ncidp, 't', varidt)
159    PRINT *, 'ncidp,varidt', ncidp, varidt
160
161    rcode = nf90_inq_varid(ncidp, 'mfu', varidmfu)
162    PRINT *, 'ncidp,varidmfu', ncidp, varidmfu
163
164    rcode = nf90_inq_varid(ncidp, 'mfd', varidmfd)
165    PRINT *, 'ncidp,varidmfd', ncidp, varidmfd
166
167    rcode = nf90_inq_varid(ncidp, 'en_u', varidenu)
168    PRINT *, 'ncidp,varidenu', ncidp, varidenu
169
170    rcode = nf90_inq_varid(ncidp, 'de_u', variddeu)
171    PRINT *, 'ncidp,variddeu', ncidp, variddeu
172
173    rcode = nf90_inq_varid(ncidp, 'en_d', varidend)
174    PRINT *, 'ncidp,varidend', ncidp, varidend
175
176    rcode = nf90_inq_varid(ncidp, 'de_d', varidded)
177    PRINT *, 'ncidp,varidded', ncidp, varidded
178
179    rcode = nf90_inq_varid(ncidp, 'coefh', varidch)
180    PRINT *, 'ncidp,varidch', ncidp, varidch
181
182    ! Thermiques
183    rcode = nf90_inq_varid(ncidp, 'fm_th', varidfmth)
184    PRINT *, 'ncidp,varidfmth', ncidp, varidfmth
185
186    rcode = nf90_inq_varid(ncidp, 'en_th', varidenth)
187    PRINT *, 'ncidp,varidenth', ncidp, varidenth
188
189    rcode = nf90_inq_varid(ncidp, 'frac_impa', varidfi)
190    PRINT *, 'ncidp,varidfi', ncidp, varidfi
191
192    rcode = nf90_inq_varid(ncidp, 'frac_nucl', varidfn)
193    PRINT *, 'ncidp,varidfn', ncidp, varidfn
194
195    rcode = nf90_inq_varid(ncidp, 'pyu1', varidyu1)
196    PRINT *, 'ncidp,varidyu1', ncidp, varidyu1
197
198    rcode = nf90_inq_varid(ncidp, 'pyv1', varidyv1)
199    PRINT *, 'ncidp,varidyv1', ncidp, varidyv1
200
201    rcode = nf90_inq_varid(ncidp, 'ftsol1', varidfts1)
202    PRINT *, 'ncidp,varidfts1', ncidp, varidfts1
203
204    rcode = nf90_inq_varid(ncidp, 'ftsol2', varidfts2)
205    PRINT *, 'ncidp,varidfts2', ncidp, varidfts2
206
207    rcode = nf90_inq_varid(ncidp, 'ftsol3', varidfts3)
208    PRINT *, 'ncidp,varidfts3', ncidp, varidfts3
209
210    rcode = nf90_inq_varid(ncidp, 'ftsol4', varidfts4)
211    PRINT *, 'ncidp,varidfts4', ncidp, varidfts4
212
213    rcode = nf90_inq_varid(ncidp, 'psrf1', varidpsr1)
214    PRINT *, 'ncidp,varidpsr1', ncidp, varidpsr1
215
216    rcode = nf90_inq_varid(ncidp, 'psrf2', varidpsr2)
217    PRINT *, 'ncidp,varidpsr2', ncidp, varidpsr2
218
219    rcode = nf90_inq_varid(ncidp, 'psrf3', varidpsr3)
220    PRINT *, 'ncidp,varidpsr3', ncidp, varidpsr3
221
222    rcode = nf90_inq_varid(ncidp, 'psrf4', varidpsr4)
223    PRINT *, 'ncidp,varidpsr4', ncidp, varidpsr4
224
225    ! ID pour les dimensions
226
227    status = nf90_inq_dimid(ncidp, 'y', yid)
228    status = nf90_inq_dimid(ncidp, 'x', xid)
229    status = nf90_inq_dimid(ncidp, 'sig_s', zid)
230    status = nf90_inq_dimid(ncidp, 'time_counter', tid)
231
232    ! lecture des dimensions
233
234    status = nf90_inquire_dimension(ncidp, yid, namedim, ncjm)
235    status = nf90_inquire_dimension(ncidp, xid, namedim, ncim)
236    status = nf90_inquire_dimension(ncidp, zid, namedim, nckev)
237    status = nf90_inquire_dimension(ncidp, tid, namedim, ncrec)
238
239    zrec = ncrec
240    zkev = nckev
241    zim = ncim
242    zjm = ncjm
243
244    zkon = zim*(zjm-2) + 2
245
246    WRITE (*, *) 'read_pstoke : zrec = ', zrec
247    WRITE (*, *) 'read_pstoke : kev = ', zkev
248    WRITE (*, *) 'read_pstoke : zim = ', zim
249    WRITE (*, *) 'read_pstoke : zjm = ', zjm
250    WRITE (*, *) 'read_pstoke : kon = ', zkon
251
252    ! niveaux de pression
253
254    status = nf90_get_var(ncidp, varidpl, pl, [1], [kev])
255
256    ! lecture de aire et phis
257
258    start(1) = 1
259    start(2) = 1
260    start(3) = 1
261    start(4) = 0
262
263    count(1) = zim
264    count(2) = zjm
265    count(3) = 1
266    count(4) = 0
267
268
269    ! **** Geopotentiel au sol ***************************************
270    ! phis
271    status = nf90_get_var(ncidp, varidps, phisfi2, start, count)
272    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, phisfi2, phisfi)
273
274    ! **** Aires des mails aux sol ************************************
275    ! aire
276    status = nf90_get_var(ncidp, varidai, airefi2, start, count)
277    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, airefi2, airefi)
278  ELSE
279
280    PRINT *, 'ok1'
281
282    ! ---------------------
283    ! lecture des champs
284    ! ---------------------
285
286    PRINT *, 'WARNING!!! Il n y a pas de test de coherence'
287    PRINT *, 'sur le nombre de niveaux verticaux dans le fichier nc'
288
289    start(1) = 1
290    start(2) = 1
291    start(3) = 1
292    start(4) = irec
293
294    count(1) = zim
295    count(2) = zjm
296    count(3) = kev
297    count(4) = 1
298
299    ! **** Temperature ********************************************
300    ! A FAIRE : Es-ce necessaire ?
301
302    ! abder t
303    status = nf90_get_var(ncidp, varidt, t2, start, count)
304    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, t2, t)
305
306    ! **** Flux pour la convection (Tiedtk)
307    ! ********************************************
308    ! mfu
309    status = nf90_get_var(ncidp, varidmfu, mfu2, start, count)
310    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, mfu2, mfu)
311
312    ! mfd
313    status = nf90_get_var(ncidp, varidmfd, mfd2, start, count)
314    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, mfd2, mfd)
315
316    ! en_u
317    status = nf90_get_var(ncidp, varidenu, en_u2, start, count)
318    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, en_u2, en_u)
319
320    ! de_u
321    status = nf90_get_var(ncidp, variddeu, de_u2, start, count)
322    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, de_u2, de_u)
323
324    ! en_d
325    status = nf90_get_var(ncidp, varidend, en_d2, start, count)
326    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, en_d2, en_d)
327
328    ! de_d
329    status = nf90_get_var(ncidp, varidded, de_d2, start, count)
330    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, de_d2, de_d)
331
332    ! **** Coefficient de mellange turbulent
333    ! *******************************************
334    ! coefh
335    PRINT *, 'LECTURE de coefh a irec =', irec
336    status = nf90_get_var(ncidp, varidch, coefh2, start, count)
337    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, coefh2, coefh)
338    ! CALL dump2d(iip1,jjp1,coefh2(1,2),'COEFH2READ   ')
339    ! CALL dump2d(iim ,jjm ,coefh (2,2),'COEFH2READ   ')
340
341    ! **** Flux ascendants et entrant dans le thermique
342    ! **********************************
343    ! Thermiques
344    PRINT *, 'LECTURE de fm_therm a irec =', irec
345    status = nf90_get_var(ncidp, varidfmth, fm_therm2, start, count)
346    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, fm_therm2, fm_therm)
347    PRINT *, 'LECTURE de en_therm a irec =', irec
348    status = nf90_get_var(ncidp, varidenth, en_therm2, start, count)
349    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, en_therm2, en_therm)
350
351    ! **** Coefficients de lessivage
352    ! *******************************************
353    ! frac_impa
354    status = nf90_get_var(ncidp, varidfi, frac_impa2, start, count)
355    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, frac_impa2, frac_impa)
356
357    ! frac_nucl
358
359    status = nf90_get_var(ncidp, varidfn, frac_nucl2, start, count)
360    CALL gr_ecrit_fi(kev, kon, nbp_lon, nbp_lat, frac_nucl2, frac_nucl)
361
362    ! **** Vents aux sol ********************************************
363
364    start(3) = irec
365    start(4) = 0
366    count(3) = 1
367    count(4) = 0
368
369    ! pyu1
370    PRINT *, 'LECTURE de yu1 a irec =', irec
371    status = nf90_get_var(ncidp, varidyu1, pyu12, start, count)
372    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, pyu12, pyu1)
373
374    ! pyv1
375    PRINT *, 'LECTURE de yv1 a irec =', irec
376    status = nf90_get_var(ncidp, varidyv1, pyv12, start, count)
377    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, pyv12, pyv1)
378
379    ! **** Temerature au sol ********************************************
380    ! ftsol1
381    PRINT *, 'LECTURE de ftsol1 a irec =', irec
382    status = nf90_get_var(ncidp, varidfts1, ftsol12, start, count)
383    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, ftsol12, ftsol1)
384
385    ! ftsol2
386    PRINT *, 'LECTURE de ftsol2 a irec =', irec
387    status = nf90_get_var(ncidp, varidfts2, ftsol22, start, count)
388    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, ftsol22, ftsol2)
389
390    ! ftsol3
391    PRINT *, 'LECTURE de ftsol3 a irec =', irec
392    status = nf90_get_var(ncidp, varidfts3, ftsol32, start, count)
393    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, ftsol32, ftsol3)
394
395    ! ftsol4
396    status = nf90_get_var(ncidp, varidfts4, ftsol42, start, count)
397    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, ftsol42, ftsol4)
398
399    ! **** Nature sol ********************************************
400    ! psrf1
401    status = nf90_get_var(ncidp, varidpsr1, psrf12, start, count)
402    ! CALL dump2d(iip1-1,jjm+1,psrf12,'PSRF1NC')
403    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, psrf12, psrf1)
404
405    ! psrf2
406    status = nf90_get_var(ncidp, varidpsr2, psrf22, start, count)
407    ! CALL dump2d(iip1-1,jjm+1,psrf22,'PSRF2NC')
408    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, psrf22, psrf2)
409
410    ! psrf3
411    status = nf90_get_var(ncidp, varidpsr3, psrf32, start, count)
412    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, psrf32, psrf3)
413
414    ! psrf4
415    status = nf90_get_var(ncidp, varidpsr4, psrf42, start, count)
416    CALL gr_ecrit_fi(1, kon, nbp_lon, nbp_lat, psrf42, psrf4)
417
418    DO i = 1, kon
419
420      psrf(i, 1) = psrf1(i)
421      psrf(i, 2) = psrf2(i)
422      psrf(i, 3) = psrf3(i)
423      ! test abderr
424      ! PRINT*,'Dans read_pstoke psrf3 =',psrf3(i),i
425      psrf(i, 4) = psrf4(i)
426
427      ftsol(i, 1) = ftsol1(i)
428      ftsol(i, 2) = ftsol2(i)
429      ftsol(i, 3) = ftsol3(i)
430      ftsol(i, 4) = ftsol4(i)
431
432    END DO
433
434  END IF
435
436  RETURN
437
438END SUBROUTINE read_pstoke0
439
Note: See TracBrowser for help on using the repository browser.