ETS  \$Id: Doxyfile 2162 2020-02-26 14:16:09Z g2dpc $
 All Classes Files Functions Variables Pages
cf_ini.f
Go to the documentation of this file.
1  SUBROUTINE cf_init( k_auto, nstep, dt, time,
2  * voltpf, d_pf_mat,d_tcam_mat )
3 
4  use durs_d_modul
5 
6  include 'double.inc'
7  !include 'param.inc'
8  include 'prm.inc'
9  include 'comevl.inc'
10  include 'parcur.inc'
11 
12 c -----------------------------
13 
14  common/comeqg/ ncequi
15 
16 ! common /comsta/ epsro,platok,
17 ! * eqdfn,betplx,tokf,psax,b0,r0,
18 ! * alf0,alf1,alf2,bet0,bet1,bet2,rax,zax,
19 ! * nurs,igdf,n_psi,n_tht,i_bsh,keyctr,i_betp,i_eqdsk
20 
21  common /comsta/ platok,eqdfn,i_bsh
22 
23  common /c_kpr/ kpr
24  common/com_flag/ kastr
25  common/com_xwx/ kxwx
26  common/com_snf/ksnf
27 c-----------------------------------------------------------------------
28 
29  dimension rc(nclim), zc(nclim), pc(nclim), psip(nclim),
30  * vc(nclim), hc(nclim)
31 
32  dimension rc1(nclim), zc1(nclim), rc2(nclim), zc2(nclim),
33  * rc3(nclim), zc3(nclim), rc4(nclim), zc4(nclim)
34 
35  INTEGER necon(nilim), ntype(nclim)
36  dimension wecon(nilim)
37 
38  common/comst0/ res(njlim), volk(njlim), volkp1(njlim)
39 
40  common/comst1/ pjk(njlim),pjkp1(njlim),pjkp(njlim),pjkd(njlim)
41  common/comst2/ psk(njlim), pskp1(njlim),pskp(njlim),pskm1(njlim)
42  include 'parloo.inc'
43  common/comloo/ rloop(nloopp),zloop(nloopp),
44  & rprob(nprobp),zprob(nprobp),fiprob(nprobp),nloop,nprob
45 
46  dimension ccurx(nclim),ccury(nclim)
47 c***********************************************************************
48 
49  real*8 voltpf(*)
50  real*8 d_pf_mat(*)
51  real*8 d_tcam_mat(*)
52 
53 c***********************************************************************
54 
55  character*40 eqdfn
56  1149 format(a40)
57  kstep = nstep
58 
59  ninf=1
60 
61  kxwx=2 ! = 0 - without prescrib. x-points;
62  ! = 1 - with one prescrib. x-point;
63  ! = 2 - with two prescrib. x-points;
64 
65  ksnf=0 ! if ksnf=1,then kxwx must be =1
66 
67 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
68  if(k_auto.eq.0) go to 2005
69 
70 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
71  if(kastr.eq.1) go to 1150
72 
73  write(fname,'(a,a)') path(1:kname),'durs.dat'
74  open(1,file=fname,form='formatted')
75  !open(1,file='durs.dat')
76 
77  read(1,*) n_tht
78  read(1,*) n_psi
79  read(1,*) igdf
80  read(1,*) epsro
81 
82  read(1,*) nurs
83  read(1,*) keyctr
84 
85  read(1,*) i_eqdsk
86  read(1,1149) eqdfn
87 
88  read(1,*) i_betp
89  read(1,*) betplx
90  read(1,*) tokf
91  read(1,*) psax
92  read(1,*) b0,r0
93  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
94  read(1,*) alf0 !!
95  read(1,*) alf1 !! P'=alf0*(1-(1-psi)**alf1)**alf2
96  read(1,*) alf2 !!
97  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
98  read(1,*) bet0 !!
99  read(1,*) bet1 !! FF'=bet0*(1-(1-psi)**bet1)**bet2
100  read(1,*) bet2 !!
101  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
102  read(1,*) alw0 !!
103  read(1,*) alw1 !! (ro*w**2)'=alw0*(1-(1-psi)**alw1)**alw2
104  read(1,*) alw2 !!
105  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
106  read(1,*) rax
107  read(1,*) zax
108 
109  close(1)
110  if(i_eqdsk.eq.1) then
111 
112  call tab_efit(tokf,psax,eqdfn,rax,zax,b0,r0)
113 
114  nurs = -3999
115  i_betp = 0
116 
117  endif
118 
119 ! write(fname,'(a,a)') path(1:kname),'durs_d.dat'
120 ! open(1,file=fname,form='formatted')
121 ! !open(1,file='durs_d.dat')
122 ! write(1,*) n_tht,n_psi,igdf,nurs,keyctr,i_eqdsk,i_betp
123 ! write(1,*) epsro,betplx,tokf,psax,b0,r0,rax,zax
124 ! write(1,*) alf0,alf1,alf2,bet0,bet1,bet2
125 ! close(1)
126 
127 
128  1150 continue
129 
130 c***********************************************************************
131 c--- input of positions of "PF_PROBE" points:
132 c
133  CALL propnt( nout, nter, ninfw, ngra1,
134  * nprob, rprob, zprob, fiprob )
135 c***********************************************************************
136 c--- input of positions of "FL_LOOP" points:
137 c
138  CALL loopnt( nout, nter, ninfw, ngra1,
139  * nloop, rloop, zloop )
140 c***********************************************************************
141 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
142 
143 c***********************************************************************
144 c--- input parameters of pfc system and passiv conductors
145 c
146  CALL conduc( nc, ncequi, ncpfc, nfw, nbp, nvv,
147  * rc,zc, pc, vc,hc, ntype,
148  * rc1,zc1, rc2,zc2, rc3,zc3, rc4,zc4,
149  * res, volk, volkp1,
150  * necon, wecon,
151  * nout, nter, ninfw, ngra1 )
152 
153 
154 
155 c
156 c***********************************************************************
157 c--- definition induct. and selfinduct. matrix
158 c--- for "EDDY" conductors: "PPIND" from COMMON /ppidps/
159 c *******
160 c
161  CALL l_matr( nout, nter, nc, ncpfc,
162  * ntype, rc, zc, vc, hc,
163  * necon,wecon )
164 c
165  write(*,*) '**'
166  DO l=1,ncequi
167  IF( l.LE.nequi ) THEN
168  pjk(l) = pfceqw(l)
169  ELSE
170  pjk(l) = pc(ncpfc+l-nequi)
171  END IF
172  enddo
173 
174  write(fname,'(a,a)') path(1:kname),'currents.wr'
175  open(1,file=fname,form='formatted')
176  !open(1,file='currents.wr',form='formatted')
177  write(1,*) nequi,ncequi
178  write(1,*)(pjk(j),j=1,ncequi)
179  close(1)
180 
181  write(fname,'(a,a)') path(1:kname),'pfcurr.wr'
182  open(1,file=fname,form='formatted')
183  !open(1,file='currents.wr',form='formatted')
184  write(1,*) nequi
185  write(1,*)(pfceqw(j),j=1,nequi)
186  close(1)
187 
188  write(fname,'(a,a)') path(1:kname),'res_mat.wr'
189  open(1,file=fname,form='formatted')
190  !open(1,file='res_mat.wr')
191  write(1,*) (res(j),j=1,ncequi)
192  close(1)
193 
194 c
195 c..... toroidal currents of passive conductor structures
196 c
197 ! CALL PASCUR( NOUT, NTER, NEQUI, NFW, NBP, NVV,
198 ! * PJK, fwcurr, bpcurr, vvcurr )
199 c
200 c********************************************************************
201 
202  call wrcoil(nc,ncpfc,rc,zc,pc,necon,wecon)
203 
204  ngrid=1
205 
206  call auto(rc,zc,pc,nc,nstep,ngrid,
207  * rc1,zc1, rc2,zc2,
208  * rc3,zc3, rc4,zc4,
209  * ntype, necon, wecon )
210 
211  2005 continue
212 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
213  if(kastr.eq.0) then
214  write(fname,'(a,a)') path(1:kname),'inpol.dat'
215  open(1,file=fname,form='formatted')
216  !open(1,file='inpol.dat')
217  read(1,*) i_bsh
218  close(1)
219  else
220  i_bsh=1
221  endif
222 
223 ! write(fname,'(a,a)') path(1:kname),'durs_d.dat'
224 ! open(1,file=fname,form='formatted')
225 ! !open(1,file='durs_d.dat')
226 ! read(1,*) n_tht,n_psi,igdf,nurs,keyctr,i_eqdsk,i_betp
227 ! read(1,*) epsro,betplx,tokf,psax,b0,r0,rax,zax
228 ! read(1,*) alf0,alf1,alf2,bet0,bet1,bet2
229 ! close(1)
230 
231 
232 !!:::::::::::::::::__________________________________:::::::::::::::
233 
234  call noauto
235 
236  call rd_ppind
237 
238  call rd_prob( nprob, rprob, zprob, fiprob )
239 
240  call rd_loop( nloop, rloop, zloop )
241 
242 
243  write(fname,'(a,a)') path(1:kname),'currents.wr'
244  open(1,file=fname,form='formatted')
245  !open(1,file='currents.wr',form='formatted')
246  read(1,*) nequi,ncequi
247  read(1,*)(pjk(j),j=1,ncequi)
248  close(1)
249 
250  write(fname,'(a,a)') path(1:kname),'pfcurr.wr'
251  open(1,file=fname,form='formatted')
252  !open(1,file='currents.wr',form='formatted')
253  read(1,*) nequi
254  read(1,*)(pfceqw(j),j=1,nequi)
255  close(1)
256 
257  write(fname,'(a,a)') path(1:kname),'res_mat.wr'
258  open(1,file=fname,form='formatted')
259  !open(1,file='res_mat.wr')
260  read(1,*) (res(j),j=1,ncequi)
261  close(1)
262 
263  call rdcoil(nc,ncpfc,rc,zc,pc,necon,wecon)
264 
265 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
266 
267  platok = tokf
268  nstep = 0
269 
270 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
271 
272  call eqb( alf0,alf1,alf2, bet0,bet1,bet2, alw0,alw1,alw2,
273  * betplx, i_betp,
274  * keyctr, nstep, platok, rax,zax, b0,r0, psax, igdf,
275  * n_tht, n_psi, epsro, nurs, i_eqdsk,i_bsh,
276  * psi_bnd,psi0_bnd)
277  !pause 'pause'
278 
279 !!!!!!!!!!!!!!!!!!!!!!!!!!ttt!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
280 !
281 !!!!!!!!!!!!!!!!!!!!!!!!!!ttt!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
282 
283 
284 
285 
286 ! nstep = 1
287 ! keyctr =1
288 
289 ! call eqb( alf0,alf1,alf2, bet0,bet1,bet2, alw0,alw1,alw2,
290 ! * betplx, i_betp,
291 ! * keyctr, nstep, platok, rax,zax, b0,r0, psax, igdf,
292 ! * n_tht, n_psi, epsro, nurs, i_eqdsk,i_bsh,
293 ! * psi_bnd,e_psi)
294 !
295  !psi_bnd=0.d0
296 
297  call prefit(rc,zc,ncpfc,necon,wecon,rax,zax,alp_b,psi_bnd)
298 
299  do iq=1,nequi
300  !curref(iq)=PFCEQW(iq)
301  curref(iq)=0.d0
302  enddo
303 
304  call curfit_(rc,zc,ncpfc,necon,wecon,psi_bnd)
305 
306  do ik=1,nequi
307  ccurx(ik)=pfceqw(ik)
308  pjk(ik)=pfceqw(ik)
309  enddo
310 
311 
312  ! pause 'pause'
313 
314 !!!!>>>>>>symmetr.case
315  isymm=0
316 !!!!>>>>>>symmetr.case
317 
318  nursb=nurs
319  psicen=psax
320  ngav1=0
321  ftok=tokf
322  nstep=0
323  ngrid=1
324 
325 
326 ! call auto(rk,zk,tk,nk,nstep,ngrid,
327 ! * rk1,zk1, rk2,zk2,
328 ! * rk3,zk3, rk4,zk4,
329 ! * ntipe, necon, wecon )
330 !
331 
332  call eq_0(pjk,psk,ncequi,nstep,ngrid,
333  * alf0, alf1, alf2, bet0, bet1, bet2,
334  * betpol, betplx, zli3,
335  * ngav1,
336  * ftok, tokout, psicen, pscout,
337  * nursb,psi_bnd,alp_b,rax,zax,n_ctrl,b0,r0 )
338 
339  call rdexf(ncequi)
340 
341  call eq_ax( pjk, psk, ncequi, nstep,ngrid,
342  * alf0, alf1, alf2, bet0, bet1, bet2,
343  * betpol, betplx, zli3,
344  * ngav1,
345  * ftok,tokout,psicen,pscout,
346  * ereve0, erps,
347  * psi_bnd,alp_b,rax,zax,isymm )
348 
349  call wrd !
350 
351  call prefit(rc,zc,ncpfc,necon,wecon,rax,zax,alp_b,psi_bnd)
352 
353  ngrid=1
354  ereve0=1.d-7
355 
356 
357  415 nstep=nstep+1
358 
359 !!!temporary!!!!only for tin
360  !if(ERPS.lt.EREVE0*10) then
361  ! iq_cs=5
362  ! do iq=1,NEQUI
363  !curref(iq)=PFCEQW(iq)
364  !if(iq.ne.iq_cs) curref(iq)=0.d0
365  ! enddo
366  !endif
367 !!!temporary!!!!only for tin
368 
369  if(kxwx.eq.1)then
370  call precal(rc,zc,ncpfc,necon, wecon )
371  elseif(kxwx.eq.0)then
372  call precal_wx(rc,zc,ncpfc,necon, wecon )
373  elseif(kxwx.eq.2)then
374  call precal(rc,zc,ncpfc,necon, wecon )
375  endif
376 
377  if(n_ctrl.eq.1 .OR. kastr.eq.1) then !limiter point
378  call curfit_l(rc,zc,ncpfc,necon,wecon,psi_bnd)
379  else
380  call curfit(rc,zc,ncpfc,necon,wecon,psi_bnd)
381  !call curfit_L(rc,zc,ncpfc,NECON,WECON,psi_bnd)
382  endif
383 
384 !!!!!!!!!!!!!!!!!!!!!saturation!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
385 ! pf_sat=4.d-3
386 ! if(nstep/5*5.eq.nstep) then
387 ! do ik=1,nequi
388 ! if(dabs(PFCEQW(ik)).gt.pf_sat) then
389 ! d_wght(ik)=d_wght(ik)*2.0d0
390 ! endif
391 ! !PFCEQW(ik)=(ccurx(ik)+PFCEQW(ik))*0.5d0
392 ! enddo
393 ! endif
394 !!!!!!!!!!!!!!!!!!!!!saturation!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
395 
396  do ik=1,nequi
397  ccurx(ik)=pfceqw(ik)
398  pjk(ik)=pfceqw(ik)
399  enddo
400 
401  call eq_ax( pjk, psk, ncequi, nstep,ngrid,
402  * alf0, alf1, alf2, bet0, bet1, bet2,
403  * betpol, betplx, zli3,
404  * ngav1,
405  * ftok,tokout,psicen,pscout,
406  * ereve0, erps,
407  * psi_bnd,alp_b,rax,zax,isymm )
408 
409  call wrd
410 
411  if(erps.gt.ereve0) go to 415
412  call wrd
413 
414  write(fname,'(a,a)') path(1:kname),'currents.wr'
415  open(1,file=fname,form='formatted')
416  !open(1,file='currents.wr',form='formatted')
417  write(1,*) nequi,ncequi
418  write(1,*)(pjk(j),j=1,ncequi)
419  close(1)
420 
421  write(fname,'(a,a)') path(1:kname),'pfc_curr.wr'
422  open(1,file=fname,form='formatted')
423  write(1,*) 'coil currents [mA]'
424  write(1,'(a6,e13.5)') 'PF1',pjk(1) !*1.d3
425  write(1,'(a6,e13.5)') 'PF2',pjk(2) !*1.d3
426  write(1,'(a6,e13.5)') 'F3,4',pjk(3) !*1.d3
427  write(1,'(a6,e13.5)') 'F5,6',pjk(4) !*1.d3
428  write(1,'(a6,e13.5)') 'CS',pjk(5) !*1.d3
429  !write(1,'(a6,e13.5)') 'CS1',pjk(1) !*1.d3
430  !write(1,'(a6,e13.5)') 'CS2U',pjk(2) !*1.d3
431  !write(1,'(a6,e13.5)') 'CS2L',pjk(3) !*1.d3
432  !write(1,'(a6,e13.5)') 'CS3U',pjk(4) !*1.d3
433  !write(1,'(a6,e13.5)') 'CS3L',pjk(5) !*1.d3
434  !write(1,'(a6,e13.5)') 'PF1',pjk(6) !*1.d3
435  !write(1,'(a6,e13.5)') 'PF2',pjk(7) !*1.d3
436  !write(1,'(a6,e13.5)') 'PF3',pjk(8) !*1.d3
437  !write(1,'(a6,e13.5)') 'PF4',pjk(9) !*1.d3
438  !write(1,'(a6,e13.5)') 'PF5',pjk(10) !*1.d3
439  !write(1,'(a6,e13.5)') 'PF6',pjk(11) !*1.d3
440  !write(1,'(a6,e13.5)') ' ',pjk(12) !*1.d3
441  !write(1,'(a6,e13.5)') ' ',pjk(13) !*1.d3
442  !write(1,'(a6,e13.5)') ' ',pjk(14) !*1.d3
443  !write(1,'(a6,e13.5)') ' ',pjk(15) !*1.d3
444  !write(1,'(a6,e13.5)') ' ',pjk(16) !*1.d3
445  !write(1,'(a6,e13.5)') ' ',pjk(17) !*1.d3
446  !write(1,'(a6,e13.5)') ' ',pjk(18) !*1.d3
447  close(1)
448 
449  write(fname,'(a,a)') path(1:kname),'tcurrs.wr'
450  open(1,file=fname,form='formatted')
451  !open(1,file='tcurrs.wr')
452  write(*,*) 'coil currents'
453 
454  do ik=1,nequi
455  do j=1,npfc
456  if(ik .eq. nepfc(j)) then
457  write(*,'(i4,2(1pe13.5))') ik,pfcw1(j),pfceqw(ik)
458  exit
459  endif
460  enddo
461  enddo
462  !call output(ngrid,betpol,zli3)
463  !call wrdcur
464  !call accr_ch
465 c write(nout,'(8H li(3)=,1pe13.5)') zli3
466 
467  !close(NOUT)
468  close(1)
469 
470  !!! stop !!!!!!!
471  pause 'pause'
472 
473  i_bsh=-1
474  nstep = 0
475  keyctr =0
476  i_betp =0
477  i_eqdsk =0
478 
479  call eqb( alf0,alf1,alf2, bet0,bet1,bet2, alw0,alw1,alw2,
480  * betplx, i_betp,
481  * keyctr, nstep, platok, rax,zax, b0,r0, psax, igdf,
482  * n_tht, n_psi, epsro, nurs, i_eqdsk,i_bsh,
483  * psi_bnd,e_psi)
484 
485 
486  !call wr_spik
487 
488  return !!!!!!!!!!
489  end
490 
491 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
492 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
493  subroutine aspid_flag(k_astr)
494  common/com_flag/kastr
495  kastr=k_astr
496  return
497  end
498 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
499  subroutine aspid_ini(k_ini)
500  common/com_ini/key_ini
501  key_ini=k_ini
502  return
503  end
tab_efit
subroutine tab_efit(tokf, psax, eqdfn, rax, zax, b0, r0)
Definition: B_wrd.f:447
conduc
subroutine conduc(NC, NCEQUI, NCPFC, NFW, NBP, NVV, RC, ZC, PC, VC, HC, NTYPE, RC1, ZC1, RC2, ZC2, RC3, ZC3, RC4, ZC4, RES, VOLK, VOLKP1, NECON, WECON, NOUT, NTER, NINFW, ngra1)
Definition: Ev1.f:826
loopnt
subroutine loopnt(NOUT, NTER, NINFW, NGRA1,
Definition: Ev3.f:254
rd_ppind
subroutine rd_ppind
Definition: Ev3.f:166
rdcoil
subroutine rdcoil(nk, nkcoil, rk, zk, tk, necon, wecon)
Definition: wrc.f:23
eqb
subroutine eqb(alf0, alf1, alf2, bet0, bet1, bet2, alw0, alw1, alw2,
Definition: B_eqb.f:3
rd_loop
subroutine rd_loop(NLOO, RLOO, ZLOO)
Definition: Ev3.f:302
l_matr
subroutine l_matr(NOUT,NTER, NC, NCPFC,
Definition: Ev3.f:63
noauto
subroutine noauto
Definition: Eq_m.f:105
auto
subroutine auto(rk, zk, tk, nk, nstep, ngrid,
Definition: Eq_m.f:2
wrd
subroutine wrd
Definition: Eq2_m.f:1378
eq_ax
subroutine eq_ax(pcequi, psicon, ncequi, nstep, ngrid,
Definition: Eq_m.f:756
wrcoil
subroutine wrcoil(nk, nkcoil, rk, zk, tk, necon, wecon)
Definition: wrc.f:1
curfit
subroutine curfit(iopt, m, x, y, w, xb, xe, k, s, nest, n, t, c, fp, wrk, lwrk, iwrk, ier)
Definition: curfit.f:1
eq_0
subroutine eq_0(pcequi, psitok, ncequi, nstep, ngrid,
Definition: Eq_m.f:148
propnt
subroutine propnt(NOUT, NTER, NINFW, NGRA1,
Definition: Ev3.f:187
rd_prob
subroutine rd_prob(NPRO, RPRO, ZPRO, FIPRO)
Definition: Ev3.f:233
rdexf
subroutine rdexf(ncequi)
Definition: EQ1_m.f:509
matrix
subroutine matrix
Definition: Matrix.f:2