ETS  \$Id: Doxyfile 2162 2020-02-26 14:16:09Z g2dpc $
 All Classes Files Functions Variables Pages
gausian_src.f90
Go to the documentation of this file.
1 MODULE gausian_src
2 
3 
4 CONTAINS
5 !-------------------------------------------------------!
6 !-------------------------------------------------------!
7 
8  SUBROUTINE gausian_sources (COREPROF, EQUILIBRIUM, CORESOURCE, code_parameters)
9 
10 !-------------------------------------------------------!
11 ! This routine provides dummy source for the !
12 ! ETS workflow. !
13 !-------------------------------------------------------!
14 ! Source: --- !
15 ! Developers: D.Kalupin !
16 ! Kontacts: Denis.Kalupin@efda.org !
17 ! !
18 ! Comments: input parameter list is specified !
19 ! in "source_dummy.xml" file. !
20 ! !
21 ! output CORESOURCE CPO is !
22 ! allocated inside the module !
23 ! !
24 !-------------------------------------------------------!
25 
26 
27  USE euitm_schemas
28  USE euitm_routines
29  USE allocate_deallocate
30  USE itm_types
31  USE copy_structures
32  USE deallocate_structures
33 
34  IMPLICIT NONE
35 
36 
37  INTEGER :: ifail
38 
39 
40 ! +++ CPO derived types:
41  TYPE (type_equilibrium), POINTER :: equilibrium(:) !input CPO with geometry quantities from previous time
42  TYPE (type_coreprof), POINTER :: coreprof(:) !input CPO with internal ETS parameters profiles from previous time
43  TYPE (type_coresource), POINTER :: coresource(:) !output CPO with sources
44  TYPE (type_param) :: code_parameters
45 
46 
47 ! +++ Local variables:
48  REAL(R8) :: time
49  REAL(R8), ALLOCATABLE :: amn(:)
50  REAL(R8), ALLOCATABLE :: zn(:)
51  REAL(R8), ALLOCATABLE :: zion(:)
52  REAL(R8), ALLOCATABLE :: rho(:)
53  REAL(R8), ALLOCATABLE :: vprime(:)
54  REAL(R8), ALLOCATABLE :: jni(:)
55  REAL(R8), ALLOCATABLE :: qel(:)
56  REAL(R8), ALLOCATABLE :: qion(:,:)
57  REAL(R8), ALLOCATABLE :: sel(:)
58  REAL(R8), ALLOCATABLE :: sion(:,:)
59  REAL(R8), ALLOCATABLE :: uion(:,:)
60  REAL(R8) :: r0
61 
62 
63 ! +++ Dimensions:
64  INTEGER, PARAMETER :: nslice = 1
65  INTEGER :: npsi
66  INTEGER :: nrho
67  INTEGER :: nnucl
68  INTEGER :: nion
69  INTEGER :: nimp
70  INTEGER, ALLOCATABLE :: nzimp(:)
71  INTEGER :: nneut
72  INTEGER, ALLOCATABLE :: ncomp(:)
73  INTEGER, ALLOCATABLE :: ntype(:)
74 
75  INTEGER :: nion2
76  INTEGER :: inucl
77  INTEGER :: iion1
78  INTEGER :: iion2
79 
80 
81 ! +++ Set dimensions:
82  nrho = SIZE (coreprof(1)%rho_tor, dim=1)
83  npsi = SIZE (equilibrium(1)%profiles_1d%rho_tor, dim=1)
84  nion2 = 100
85  CALL get_comp_dimensions(coreprof(1)%COMPOSITIONS, nnucl, nion, nimp, nzimp, nneut, ntype, ncomp)
86 
87 ! +++ Allocate output CPO:
88  CALL allocate_coresource_cpo(nslice, nrho, nnucl, nion, nimp, nzimp, nneut, ntype, ncomp, coresource)
89  CALL deallocate_cpo(coresource(1)%COMPOSITIONS)
90  CALL copy_cpo(coreprof(1)%COMPOSITIONS, coresource(1)%COMPOSITIONS)
91 
92 
93 ! +++ Allocate local variables:
94  ALLOCATE ( amn(nion2) )
95  ALLOCATE ( zn(nion2) )
96  ALLOCATE ( zion(nion2) )
97  ALLOCATE ( rho(nrho) )
98  ALLOCATE ( vprime(nrho) )
99  ALLOCATE ( jni(nrho) )
100  ALLOCATE ( qel(nrho) )
101  ALLOCATE ( qion(nrho,nion2))
102  ALLOCATE ( sel(nrho) )
103  ALLOCATE ( sion(nrho,nion2))
104  ALLOCATE ( uion(nrho,nion2))
105 
106 
107  amn = 0.0_r8
108  zn = 0.0_r8
109  zion = 0.0_r8
110  rho = 0.0_r8
111  vprime = 0.0_r8
112  jni = 0.0_r8
113  qel = 0.0_r8
114  qion = 0.0_r8
115  sel = 0.0_r8
116  sion = 0.0_r8
117  uion = 0.0_r8
118 
119 
120 
121 ! +++ Save output in CPO:
122  time = coreprof(1)%time !time [s]
123 
124  rho = coreprof(1)%rho_tor !rho [m]
125 
126  r0 = equilibrium(1)%global_param%toroid_field%r0
127 
128  CALL l3deriv(equilibrium(1)%profiles_1d%volume, equilibrium(1)%profiles_1d%rho_tor, npsi, &
129  vprime, rho, nrho)
130 
131 
132 ! +++ Calculate sources:
133 ! << input >> << output >>
134  CALL additional_source(amn, zn, zion, nrho, rho, r0, vprime, nion2, qel, qion, sel, sion, uion, jni, code_parameters)
135 
136 
137 
138  nion2 = min(SIZE(amn),SIZE(zn),SIZE(zion))
139 
140 
141 ! +++ Save output in CPO:
142  coresource(1)%time = time !time [s]
143 
144  coresource(1)%VALUES(1)%rho_tor = rho !rho [m]
145 
146  coresource(1)%VALUES(1)%j = jni !j_ni [A/m^2]
147 
148  coresource(1)%VALUES(1)%qe%exp = qel !Qe_exp [W/m^3]
149  coresource(1)%VALUES(1)%qe%imp = 0.0e0_r8 !Qe_imp [1/m^3/s]
150 
151  coresource(1)%VALUES(1)%se%exp = sel !Se_exp [1/m^3/s]
152  coresource(1)%VALUES(1)%se%imp = 0.0e0_r8 !Se_imp [1/s]
153 
154  DO iion1 = 1, nion
155  inucl = coresource(1)%COMPOSITIONS%IONS(iion1)%nucindex
156  DO iion2 = 1, nion2
157  IF (abs(amn(iion2) - coresource(1)%COMPOSITIONS%NUCLEI(inucl)%amn) .LE. 0.25 .AND. &
158  abs(zn(iion2) - coresource(1)%COMPOSITIONS%NUCLEI(inucl)%zn ) .LE. 0.25 .AND. &
159  abs(zion(iion2) - coresource(1)%COMPOSITIONS%IONS(iion1)%zion ) .LE. 0.25) THEN
160 
161  coresource(1)%VALUES(1)%si%exp(:,iion1) = sion(:,iion2) !Si_exp [1/m^3/s]
162  coresource(1)%VALUES(1)%si%imp(:,iion1) = 0.0e0_r8 !Si_imp [1/s]
163 
164  coresource(1)%VALUES(1)%qi%exp(:,iion1) = qion(:,iion2) !Qi_exp [W/m^3]
165  coresource(1)%VALUES(1)%qi%imp(:,iion1) = 0.0e0_r8 !Qi_imp [1/m^3/s]
166 
167  coresource(1)%VALUES(1)%ui%exp(:,iion1) = uion(:,iion2) !Ui_exp [kg/m/s^2]
168  coresource(1)%VALUES(1)%ui%imp(:,iion1) = 0.0e0_r8 !Ui_imp [kg/m^2/s]
169 
170  END IF
171  END DO
172  END DO
173 
174 
175 
176 ! +++ ADD IDENTIFIER TO OUTPUT CPO VALUES(1):
177  ALLOCATE (coresource(1)%VALUES(1)%sourceid%id(1))
178  ALLOCATE (coresource(1)%VALUES(1)%sourceid%description(1))
179  coresource(1)%VALUES(1)%sourceid%id = 'gaussian'
180  coresource(1)%VALUES(1)%sourceid%flag = 33
181  coresource(1)%VALUES(1)%sourceid%description = 'Gaussian'
182 
183 
184 
185 ! +++ Deallocate local variables:
186  IF(ALLOCATED(amn)) DEALLOCATE (amn)
187  IF(ALLOCATED(rho)) DEALLOCATE (rho)
188  IF(ALLOCATED(vprime)) DEALLOCATE (vprime)
189  IF(ALLOCATED(jni)) DEALLOCATE (jni)
190  IF(ALLOCATED(qel)) DEALLOCATE (qel)
191  IF(ALLOCATED(qion)) DEALLOCATE (qion)
192  IF(ALLOCATED(sel)) DEALLOCATE (sel)
193  IF(ALLOCATED(sion)) DEALLOCATE (sion)
194  IF(ALLOCATED(uion)) DEALLOCATE (uion)
195  IF(ALLOCATED(nzimp)) DEALLOCATE (nzimp)
196  IF(ALLOCATED(ncomp)) DEALLOCATE (ncomp)
197  IF(ALLOCATED(ntype)) DEALLOCATE (ntype)
198 
199 
200 ! + + + + + + + + + + + + + + + + + + + + + + + + + + + +
201 ! + + + + + + + + + + + + + + + + + + + + + + + + + + + +
202 
203 
204 
205  RETURN
206 
207 
208  END SUBROUTINE gausian_sources
209 
210 !-------------------------------------------------------!
211 !-------------------------------------------------------!
212 
213 
214 
215 
216 
217 
218 !-------------------------------------------------------!
219 !-------------------------------------------------------!
220 
221  SUBROUTINE additional_source (AMN, ZN, ZION, NRHO, RHO, R0, VPRIME, NION, QEL, QION, SEL, SION, UION, JNI, code_parameters)
222 
223  USE itm_constants
224  USE itm_types
225 
226  USE euitm_routines
227  USE euitm_schemas
228  USE euitm_xml_parser
229  USE xml_file_reader
230 
231  IMPLICIT NONE
232 
233 
234 ! +++ Dimensions:
235  INTEGER :: i
236  INTEGER :: nrho, irho
237  INTEGER :: nion
238 
239 
240 ! +++ Local variables:
241  REAL(R8) :: rho(nrho)
242  REAL(R8) :: vprime(nrho)
243 
244  REAL(R8) :: amn(nion)
245  REAL(R8) :: zn(nion)
246  REAL(R8) :: zion(nion)
247 
248  REAL(R8) :: jni(nrho)
249  REAL(R8) :: qel(nrho)
250  REAL(R8) :: qion(nrho,nion)
251  REAL(R8) :: sel(nrho)
252  REAL(R8) :: sion(nrho,nion)
253  REAL(R8) :: uion(nrho,nion)
254  REAL(R8) :: r0
255 
256  REAL(R8) :: gfun(nrho)
257  REAL(R8) :: a, c, w, intfun
258  REAL(R8) :: ftot
259 
260 
261 ! +++ Parameters from XML:
262 
263 ! CURRENT:
264  REAL(R8) :: jnitot !total noninductive current [A]
265  REAL(R8) :: rcurr !rho position of current profile maximum [m]
266  REAL(R8) :: fwcurr !A full width at half maximum of the current profiles [m]
267 
268 ! ELECTRON HEATING:
269  REAL(R8) :: wtot_el !total energy input [W]
270  REAL(R8) :: rheat_el !rho position of heating profile maximum [m]
271  REAL(R8) :: fwheat_el !A full width at half maximum of the heating profiles [m]
272 
273 ! ION HEATING:
274  REAL(R8) :: wtot(nion) !total energy input [W]
275  REAL(R8) :: rheat(nion) !rho position of heating profile maximum [m]
276  REAL(R8) :: fwheat(nion) !A full width at half maximum of the heating profiles [m]
277 
278 ! ELECTRONS:
279  REAL(R8) :: stot_el !total electron particle input [s^-1]
280  REAL(R8) :: rpart_el !rho position of electron particle source profile maximum [m]
281  REAL(R8) :: fwpart_el !A full width at half maximum of the electron particle source profiles [m]
282 
283 ! IONS:
284  REAL(R8) :: stot(nion) !total ion particle input [s^-1]
285  REAL(R8) :: rpart(nion) !rho position of ion particle source profile maximum [m]
286  REAL(R8) :: fwpart(nion) !A full width at half maximum of the ion particle source profiles [m]
287 
288 ! MOMENTUM:
289  REAL(R8) :: utot(nion) !Total momentum [kg*m^2*s^-1]
290  REAL(R8) :: rmom(nion) !rho position of momentum source profile maximum [m]
291  REAL(R8) :: fwmom(nion) !A full width at half maximum of the momentum source profiles [m]
292 
293 
294 ! +++ Other
295  LOGICAL, SAVE :: first = .true.
296  INTEGER :: return_status
297  TYPE (type_param) :: code_parameters
298 
299 
300 ! CALL FILL_PARAM (code_parameters, 'XML/source_dummy.xml', '', 'XML/source_dummy.xsd')
301 
302 ! choose some defaults
303  jnitot = 0.0_r8
304  rcurr = 0.0_r8
305  fwcurr = 0.5_r8 * rho(nrho)
306 
307  wtot_el = 0.0_r8
308  rheat_el = 0.0_r8
309  fwheat_el = 0.5_r8 * rho(nrho)
310 
311  wtot = 0.0_r8
312  rheat = 0.0_r8
313  fwheat = 0.5_r8 * rho(nrho)
314 
315  stot_el = 0.0_r8
316  rpart_el = 0.0_r8
317  fwpart_el = 0.5_r8 * rho(nrho)
318 
319  stot = 0.0_r8
320  rpart = 0.0_r8
321  fwpart = 0.5_r8 * rho(nrho)
322 
323  utot = 0.0_r8
324  rmom = 0.0_r8
325  fwmom = 0.5_r8 * rho(nrho)
326 
327  CALL assign_dummy_sources(code_parameters, return_status)
328 
329  IF (return_status /= 0) THEN
330  WRITE(*,*) 'ERROR: Could not assign source multipliers.'
331  END IF
332 
333 
334 
335 !-------------------------------------------------------!
336 ! +++ Ambipolarity check:
337 
338  IF (stot_el.NE.sum(stot*zion))THEN
339  WRITE(*,*)'=================================='
340  WRITE(*,*)'=================================='
341  WRITE(*,*)'== Your Gaussian sources are =='
342  WRITE(*,*)'== not ambipolar! =='
343  WRITE(*,*)'== This can break your worflow. =='
344  WRITE(*,*)'=================================='
345  WRITE(*,*)'=================================='
346  END IF
347 
348 !-------------------------------------------------------!
349 ! +++ Current sources:
350 
351  IF(jnitot .NE. 0.0_r8) THEN
352  a = rcurr
353  c = fwcurr/2.35482_r8
354  w = jnitot
355 
356  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
357  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
358  * vprime / 2.0_r8 / itm_pi / rho
359 
360  IF (rho(1).EQ.0.0_r8) &
361  gfun(1) = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
362  * exp(-a**2 / 2.0_r8 / c**2) &
363  * 2.0_r8 * itm_pi * r0
364 
365  CALL integ(nrho, rho, gfun, intfun)
366 
367  jni = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
368  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
369  * w / intfun
370  ELSE
371  jni = 0.0_r8
372  ENDIF
373 
374 !-------------------------------------------------------!
375 ! +++ Electron heating sources:
376 
377  IF(wtot_el .NE. 0.0_r8) THEN
378  a = rheat_el
379  c = fwheat_el/2.35482_r8
380  w = wtot_el
381 
382  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
383  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
384  * vprime
385 
386  CALL integ(nrho, rho, gfun, intfun)
387 
388  qel = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
389  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
390  * w / intfun
391  ELSE
392  qel = 0.0_r8
393  ENDIF
394 
395 !-------------------------------------------------------!
396 ! +++ Individual ion heating sources:
397 
398  DO i = 1,nion
399 
400  IF(wtot(i) .NE. 0.0_r8) THEN
401  a = rheat(i)
402  c = fwheat(i)/2.35482_r8
403  w = wtot(i)
404 
405  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
406  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
407  * vprime
408 
409  CALL integ(nrho, rho, gfun, intfun)
410 
411  qion(:,i)= 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
412  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
413  * w / intfun
414  ELSE
415  qion(:,i) = 0.0_r8
416  ENDIF
417 
418  END DO
419 
420 !-------------------------------------------------------!
421 ! +++ Individual electron sources:
422 
423  IF(stot_el .NE. 0.0) THEN
424  a = rpart_el
425  c = fwpart_el/2.35482_r8
426  w = stot_el
427 
428  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
429  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
430  * vprime
431 
432  CALL integ(nrho, rho, gfun, intfun)
433 
434  sel = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
435  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
436  * w / intfun
437  ELSE
438  sel = 0.0_r8
439  ENDIF
440 
441 !-------------------------------------------------------!
442 ! +++ Individual ion particle sources:
443 
444  DO i = 1,nion
445 
446  IF(stot(i) .NE. 0.0_r8) THEN
447  a = rpart(i)
448  c = fwpart(i)/2.35482_r8
449  w = stot(i)
450 
451  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
452  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
453  * vprime
454 
455  CALL integ(nrho, rho, gfun, intfun)
456 
457  sion(:,i)= 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
458  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
459  * w / intfun
460  ELSE
461  sion(:,i) = 0.0_r8
462  ENDIF
463 
464  END DO
465 
466 !-------------------------------------------------------!
467 ! +++ Individual momentum sources:
468 
469  DO i = 1,nion
470 
471  IF(rmom(i) .NE. 0.0) THEN
472  a = rmom(i)
473  c = fwmom(i)/2.35482_r8
474  w = utot(i)
475 
476  gfun = 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
477  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
478  * vprime
479 
480  CALL integ(nrho, rho, gfun, intfun)
481 
482  uion(:,i)= 1.0_r8 / (2.0_r8*itm_pi)**0.5/ c &
483  * exp(-(rho-a)**2 / 2.0_r8 / c**2) &
484  * w / intfun
485  ELSE
486  uion(:,i) = 0.0_r8
487  ENDIF
488 
489  END DO
490 
491 
492  RETURN
493 
494 
495 !-------------------------------------------------------!
496 !-------------------------------------------------------!
497 
498 
499 
500  CONTAINS
501 !-------------------------------------------------------!
502 !-------------------------------------------------------!
503  SUBROUTINE assign_dummy_sources(codeparameters, return_status)
504 
505 !-------------------------------------------------------!
506 ! This subroutine calls the XML parser for !
507 ! the dummy sources and assign the !
508 ! resulting values to the corresponding variables !
509 !-------------------------------------------------------!
510 ! Source: --- !
511 ! Developers: D.Kalupin !
512 ! Kontacts: Denis.Kalupin@efda.org !
513 ! !
514 ! Comments: --- !
515 ! !
516 !-------------------------------------------------------!
517 
518  USE itm_types
519  USE euitm_schemas
520  USE euitm_xml_parser
521 
522  IMPLICIT NONE
523 
524 
525  TYPE(type_param) :: codeparameters
526  INTEGER(ITM_I4) :: return_status
527 
528  TYPE(tree) :: parameter_list
529  TYPE(element), POINTER :: temp_pointer
530  INTEGER(ITM_I4) :: i, nparm, n_values
531  CHARACTER(len = 132) :: cname
532 
533  INTEGER :: n_data
534 
535  return_status = 0
536 
537 !-- parse xml-string codeparameters%parameters
538 
539  WRITE(6,*) 'Calling euitm_xml_parse'
540  CALL euitm_xml_parse(codeparameters, nparm, parameter_list)
541  WRITE(6,*) 'Called euitm_xml_parse'
542 
543 !-- assign variables
544 
545  temp_pointer => parameter_list%first
546 
547  outer: DO
548  cname = char2str(temp_pointer%cname) ! necessary for AIX
549  SELECT CASE (cname)
550 
551 
552 !-- parameters overall
553  CASE ("parameters")
554  temp_pointer => temp_pointer%child
555  cycle
556 
557 !-----------------------------------------------------------
558 !-- Parameters for current source
559  CASE ("currents")
560  temp_pointer => temp_pointer%child
561  cycle
562 
563  CASE ("JNITOT")
564  IF (ALLOCATED(temp_pointer%cvalue)) &
565  CALL char2num(temp_pointer%cvalue, jnitot)
566 
567  CASE ("RCURR")
568  IF (ALLOCATED(temp_pointer%cvalue)) &
569  CALL char2num(temp_pointer%cvalue, rcurr)
570 
571  CASE ("FWCURR")
572  IF (ALLOCATED(temp_pointer%cvalue)) &
573  CALL char2num(temp_pointer%cvalue, fwcurr)
574 
575 
576 
577 !-----------------------------------------------------------
578 !-- Parameters for electrons
579  CASE ("electrons")
580  temp_pointer => temp_pointer%child
581  cycle
582 
583 !-- heating source
584  CASE ("heating_el")
585  temp_pointer => temp_pointer%child
586  cycle
587 
588  CASE ("WTOT_el")
589  IF (ALLOCATED(temp_pointer%cvalue)) &
590  CALL char2num(temp_pointer%cvalue, wtot_el)
591 
592  CASE ("RHEAT_el")
593  IF (ALLOCATED(temp_pointer%cvalue)) &
594  CALL char2num(temp_pointer%cvalue, rheat_el)
595 
596  CASE ("FWHEAT_el")
597  IF (ALLOCATED(temp_pointer%cvalue)) &
598  CALL char2num(temp_pointer%cvalue, fwheat_el)
599 
600 
601 
602 !-- Parameters for particle source
603  CASE ("particles_el")
604  temp_pointer => temp_pointer%child
605  cycle
606 
607  CASE ("STOT_el")
608  IF (ALLOCATED(temp_pointer%cvalue)) &
609  CALL char2num(temp_pointer%cvalue, stot_el)
610 
611  CASE ("RPART_el")
612  IF (ALLOCATED(temp_pointer%cvalue)) &
613  CALL char2num(temp_pointer%cvalue, rpart_el)
614 
615  CASE ("FWPART_el")
616  IF (ALLOCATED(temp_pointer%cvalue)) &
617  CALL char2num(temp_pointer%cvalue, fwpart_el)
618 
619 
620 
621 
622 !-----------------------------------------------------------
623 !-- Parameters for ions
624  CASE ("ions")
625  temp_pointer => temp_pointer%child
626  cycle
627 
628 !-- ion composition
629  CASE ("composition")
630  temp_pointer => temp_pointer%child
631  cycle
632 
633  CASE ("AMN")
634  IF (ALLOCATED(temp_pointer%cvalue)) &
635  CALL scan_str2real(char2str(temp_pointer%cvalue), amn, n_data)
636 
637  CASE ("ZN")
638  IF (ALLOCATED(temp_pointer%cvalue)) &
639  CALL scan_str2real(char2str(temp_pointer%cvalue), zn, n_data)
640 
641  CASE ("ZION")
642  IF (ALLOCATED(temp_pointer%cvalue)) &
643  CALL scan_str2real(char2str(temp_pointer%cvalue), zion, n_data)
644 
645 
646 !-- heating source
647  CASE ("heating")
648  temp_pointer => temp_pointer%child
649  cycle
650 
651  CASE ("WTOT")
652  IF (ALLOCATED(temp_pointer%cvalue)) &
653  CALL scan_str2real(char2str(temp_pointer%cvalue), wtot, n_data)
654 
655  CASE ("RHEAT")
656  IF (ALLOCATED(temp_pointer%cvalue)) &
657  CALL scan_str2real(char2str(temp_pointer%cvalue), rheat, n_data)
658 
659  CASE ("FWHEAT")
660  IF (ALLOCATED(temp_pointer%cvalue)) &
661  CALL scan_str2real(char2str(temp_pointer%cvalue), fwheat, n_data)
662 
663 
664 !-- Parameters for particle source
665  CASE ("particles")
666  temp_pointer => temp_pointer%child
667  cycle
668 
669  CASE ("STOT")
670  IF (ALLOCATED(temp_pointer%cvalue)) &
671  CALL scan_str2real(char2str(temp_pointer%cvalue), stot, n_data)
672 
673  CASE ("RPART")
674  IF (ALLOCATED(temp_pointer%cvalue)) &
675  CALL scan_str2real(char2str(temp_pointer%cvalue), rpart, n_data)
676 
677  CASE ("FWPART")
678  IF (ALLOCATED(temp_pointer%cvalue)) &
679  CALL scan_str2real(char2str(temp_pointer%cvalue), fwpart, n_data)
680 
681 
682 !-- Parameters for momentum source
683  CASE ("momentum")
684  temp_pointer => temp_pointer%child
685  cycle
686 
687  CASE ("UTOT")
688  IF (ALLOCATED(temp_pointer%cvalue)) &
689  CALL scan_str2real(char2str(temp_pointer%cvalue), utot, n_data)
690 
691  CASE ("RMOM")
692  IF (ALLOCATED(temp_pointer%cvalue)) &
693  CALL scan_str2real(char2str(temp_pointer%cvalue), rmom, n_data)
694 
695  CASE ("FWMOM")
696  IF (ALLOCATED(temp_pointer%cvalue)) &
697  CALL scan_str2real(char2str(temp_pointer%cvalue), fwmom, n_data)
698 
699 
700  CASE default
701  WRITE(*, *) 'ERROR: invalid parameter', cname
702  return_status = 1
703  EXIT
704  END SELECT
705 
706 
707  DO
708  IF (ASSOCIATED(temp_pointer%sibling)) THEN
709  temp_pointer => temp_pointer%sibling
710  EXIT
711  END IF
712  IF (ASSOCIATED(temp_pointer%parent, parameter_list%first )) &
713  EXIT outer
714  IF (ASSOCIATED(temp_pointer%parent)) THEN
715  temp_pointer => temp_pointer%parent
716  ELSE
717  WRITE(*, *) 'ERROR: broken list.'
718  RETURN
719  END IF
720  END DO
721  END DO outer
722 
723 
724 ! -- destroy tree
725  CALL destroy_xml_tree(parameter_list)
726 
727 
728  RETURN
729 
730 
731  END SUBROUTINE assign_dummy_sources
732 !-------------------------------------------------------!
733 !-------------------------------------------------------!
734 
735 
736  END SUBROUTINE additional_source
737 !-------------------------------------------------------!
738 !-------------------------------------------------------!
739 
740 
741 
742 
743 
744 
745 !-------------------------------------------------------!
746 !-------------------------------------------------------!
747 ! This subroutine calculates integral of a function y(x)
748 ! from X=0
749  SUBROUTINE integ(N,X,Y,INTY)
750 
751 
752  USE itm_types
753 
754 
755  IMPLICIT NONE
756 
757  INTEGER :: n ! number of radial points (input)
758  INTEGER :: i
759 
760  REAL (R8) :: x(n), & ! argument array (input)
761  y(n), & ! function array (input)
762  inty ! integral value (output)
763 
764  inty = 0.0_r8
765  DO i=2,n
766  inty = inty + (y(i)+y(i-1))*(x(i)-x(i-1))/2.e0_r8
767  END DO
768 
769  RETURN
770 
771  END SUBROUTINE integ
772 !-------------------------------------------------------!
773 !-------------------------------------------------------!
774 
775 
776 
777 
778 !-------------------------------------------------------!
779 !-------------------------------------------------------!
780  SUBROUTINE interpolate(ninput,rinput,finput,nout,rout,fout)
781 
782  ! +++ this subroutine interpolates profile Finput(rinput) to the new greed Fout(rout)
783 
784  USE itm_types
785 
786  IMPLICIT NONE
787 
788  INTEGER :: i,ii,jj, ninput,nout
789  ! Input: radii, function
790  REAL (R8) :: rinput(ninput),finput(ninput)
791  ! Output: radii, function
792  REAL (R8) :: rout(nout),fout(nout)
793  REAL (R8) :: f(2),r,h1,h2,y1,y2,y3,a,b,c,z
794 
795 
796  DO i=1,nout
797  r=rout(i)
798  IF(r.LE.rinput(2)) THEN
799  h1=rinput(2)-rinput(1)
800  h2=rinput(3)-rinput(2)
801  y1=finput(1)
802  y2=finput(2)
803  y3=finput(3)
804  a=((y1-y2)*h2+(y3-y2)*h1)/h1/h2/(h1+h2)
805  b=((y3-y2)*h1*h1-(y1-y2)*h2*h2)/h1/h2/(h1+h2)
806  c=y2
807  z=a*(r-rinput(2))**2+b*(r-rinput(2))+c
808  goto 100
809  ENDIF
810  IF(r.GT.rinput(2).AND.r.LE.rinput(ninput-1)) THEN
811  DO ii=3,ninput-1
812  IF(r.GT.rinput(ii-1).AND.r.LE.rinput(ii)) THEN
813  DO jj=1,2
814  h1=rinput(ii-2+jj)-rinput(ii-3+jj)
815  h2=rinput(ii-1+jj)-rinput(ii-2+jj)
816  y1=finput(ii-3+jj)
817  y2=finput(ii-2+jj)
818  y3=finput(ii-1+jj)
819  a=((y1-y2)*h2+(y3-y2)*h1)/h1/h2/(h1+h2)
820  b=((y3-y2)*h1*h1-(y1-y2)*h2*h2)/h1/h2/(h1+h2)
821  c=y2
822  f(jj)=a*(r-rinput(ii-2+jj))**2+b*(r-rinput(ii-2+jj))+c
823  ENDDO
824  z=(f(1)*(rinput(ii)-r)+f(2)*(r-rinput(ii-1)))/ &
825  (rinput(ii)-rinput(ii-1))
826  ENDIF
827  ENDDO
828  goto 100
829  ENDIF
830  IF(r.GE.rinput(ninput-1)) THEN
831  h1=rinput(ninput-1)-rinput(ninput-2)
832  h2=rinput(ninput)-rinput(ninput-1)
833  y1=finput(ninput-2)
834  y2=finput(ninput-1)
835  y3=finput(ninput)
836  a=((y1-y2)*h2+(y3-y2)*h1)/h1/h2/(h1+h2)
837  b=((y3-y2)*h1*h1-(y1-y2)*h2*h2)/h1/h2/(h1+h2)
838  c=y2
839  z=a*(r-rinput(ninput-1))**2+b*(r-rinput(ninput-1))+c
840  goto 100
841  ENDIF
842 100 fout(i)=z
843  ENDDO
844 
845 !DPC
846  IF(rinput(1).EQ.rout(1)) THEN
847  IF(finput(1).NE.fout(1)) THEN
848  WRITE(*,*) 'INTERPOLATION corrected for ',rinput(1),finput(1),fout(1)
849  fout(1)=finput(1)
850  ENDIF
851  ENDIF
852  RETURN
853 
854  END SUBROUTINE interpolate
855 !-------------------------------------------------------!
856 !-------------------------------------------------------!
857 
858 
859 
860 
861 !-------------------------------------------------------!
862 !-------------------------------------------------------!
863 
864 END MODULE gausian_src
865 
866 
867 
868 
869 
870 
871 
l3deriv
subroutine l3deriv(y_in, x_in, nr_in, dydx_out, x_out, nr_out)
Definition: l3interp.f90:59
gausian_src::additional_source
subroutine additional_source(AMN, ZN, ZION, NRHO, RHO, R0, VPRIME, NION, QEL, QION, SEL, SION, UION, JNI, code_parameters)
Definition: gausian_src.f90:221
allocate_deallocate::get_comp_dimensions
subroutine get_comp_dimensions(COMPOSITIONS, NNUCL, NION, NIMP, NZIMP, NNEUT, NTYPE, NCOMP)
Definition: allocate_deallocate.f90:1845
gausian_src
Definition: gausian_src.f90:1
gausian_src::integ
subroutine integ(N, X, Y, INTY)
Definition: gausian_src.f90:749
allocate_deallocate
This module contains routines for allocation/deallocation if CPOs used in ETS.
Definition: allocate_deallocate.f90:8
allocate_deallocate::allocate_coresource_cpo
subroutine allocate_coresource_cpo(NSLICE, NRHO, NNUCL, NION, NIMP, NZIMP, NNEUT, NTYPE, NCOMP, CORESOURCE)
This routine allocates CORESOURCE CPO.
Definition: allocate_deallocate.f90:460
itm_toolbox::interpolate
real(r8) function, dimension(size(x_out)) interpolate(x_in, y_in, x_out)
&quot;generic&quot; interpolation routine (only calls l3interp at the moment)
Definition: itm_toolbox.F90:29
assign_dummy_sources
subroutine assign_dummy_sources(codeparameters, return_status)
Definition: gausian_src.f90:503
gausian_src::gausian_sources
subroutine gausian_sources(COREPROF, EQUILIBRIUM, CORESOURCE, code_parameters)
Definition: gausian_src.f90:8