SFEMaNS  version 4.1 (work in progress)
Reference documentation for SFEMaNS
 All Classes Files Functions Variables Groups Pages
condlim_test_3.f90
Go to the documentation of this file.
2  use my_util
3  USE def_type_mesh
4  USE input_data
5 !!$ATTENTION
6 !!$Some subroutines have been commented to avoid warning messages when compiling executable.
7 !!$It can not be done in the module boundary_generic that expects all subroutines to be present.
8 !!$END ATTENTION
9 !!$ PUBLIC :: init_velocity_pressure
10 !!$ PUBLIC :: init_temperature
11 !!$ PUBLIC :: init_level_set
12 !!$ PUBLIC :: source_in_NS_momentum
13 !!$ PUBLIC :: source_in_temperature
14 !!$ PUBLIC :: source_in_level_set
15 !!$ PUBLIC :: vv_exact
16 !!$ PUBLIC :: imposed_velocity_by_penalty
17 !!$ PUBLIC :: pp_exact
18 !!$ PUBLIC :: temperature_exact
19 !!$ PUBLIC :: level_set_exact
20 !!$ PUBLIC :: penal_in_real_space
21 !!$ PUBLIC :: extension_velocity
22  PUBLIC :: vexact
23 !!$ PUBLIC :: H_B_quasi_static
24  PUBLIC :: hexact
25  PUBLIC :: phiexact
26  PUBLIC :: jexact_gauss
27  PUBLIC :: eexact_gauss
28  PUBLIC :: init_maxwell
29 !!$ PUBLIC :: mu_bar_in_fourier_space
30 !!$ PUBLIC :: grad_mu_bar_in_fourier_space
31 !!$ PUBLIC :: mu_in_real_space
32  PRIVATE
33  REAL (KIND=8), PRIVATE :: alpha=1.d0, beta=1.d0
34 
35 CONTAINS
36  !===============================================================================
37  ! Boundary conditions for Navier-Stokes
38  !===============================================================================
39 
40 !!$ !===Initialize velocity, pressure
41 !!$ SUBROUTINE init_velocity_pressure(mesh_f, mesh_c, time, dt, list_mode, &
42 !!$ un_m1, un, pn_m1, pn, phin_m1, phin)
43 !!$ IMPLICIT NONE
44 !!$ TYPE(mesh_type) :: mesh_f, mesh_c
45 !!$ REAL(KIND=8), INTENT(OUT):: time
46 !!$ REAL(KIND=8), INTENT(IN) :: dt
47 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
48 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: un_m1, un
49 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: pn_m1, pn, phin_m1, phin
50 !!$ INTEGER :: mode, i, j
51 !!$ REAL(KIND=8), DIMENSION(mesh_c%np) :: pn_m2
52 !!$
53 !!$ time = 0.d0
54 !!$ DO i= 1, SIZE(list_mode)
55 !!$ mode = list_mode(i)
56 !!$ DO j = 1, 6
57 !!$ !===velocity
58 !!$ un_m1(:,j,i) = vv_exact(j,mesh_f%rr,mode,time-dt)
59 !!$ un (:,j,i) = vv_exact(j,mesh_f%rr,mode,time)
60 !!$ END DO
61 !!$ DO j = 1, 2
62 !!$ !===pressure
63 !!$ pn_m2(:) = pp_exact(j,mesh_c%rr,mode,time-2*dt)
64 !!$ pn_m1 (:,j,i) = pp_exact(j,mesh_c%rr,mode,time-dt)
65 !!$ pn (:,j,i) = pp_exact(j,mesh_c%rr,mode,time)
66 !!$ phin_m1(:,j,i) = pn_m1(:,j,i) - pn_m2(:)
67 !!$ phin (:,j,i) = Pn (:,j,i) - pn_m1(:,j,i)
68 !!$ ENDDO
69 !!$ ENDDO
70 !!$ END SUBROUTINE init_velocity_pressure
71 
72 !!$ !===Initialize temperature
73 !!$ SUBROUTINE init_temperature(mesh, time, dt, list_mode, tempn_m1, tempn)
74 !!$ IMPLICIT NONE
75 !!$ TYPE(mesh_type) :: mesh
76 !!$ REAL(KIND=8), INTENT(OUT):: time
77 !!$ REAL(KIND=8), INTENT(IN) :: dt
78 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
79 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: tempn_m1, tempn
80 !!$ INTEGER :: mode, i, j
81 !!$
82 !!$ time = 0.d0
83 !!$ DO i= 1, SIZE(list_mode)
84 !!$ mode = list_mode(i)
85 !!$ DO j = 1, 2
86 !!$ tempn_m1(:,j,i) = temperature_exact(j, mesh%rr, mode, time-dt)
87 !!$ tempn (:,j,i) = temperature_exact(j, mesh%rr, mode, time)
88 !!$ ENDDO
89 !!$ ENDDO
90 !!$ END SUBROUTINE init_temperature
91 
92 !!$ !===Initialize level_set
93 !!$ SUBROUTINE init_level_set(vv_mesh, time, &
94 !!$ dt, list_mode, level_set_m1, level_set)
95 !!$ IMPLICIT NONE
96 !!$ TYPE(mesh_type) :: vv_mesh
97 !!$ REAL(KIND=8), INTENT(OUT):: time
98 !!$ REAL(KIND=8), INTENT(IN) :: dt
99 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
100 !!$ REAL(KIND=8), DIMENSION(:,:,:,:), INTENT(OUT):: level_set, level_set_m1
101 !!$ INTEGER :: mode, i, j, n
102 !!$
103 !!$ time = 0.d0
104 !!$ DO i= 1, SIZE(list_mode)
105 !!$ mode = list_mode(i)
106 !!$ DO j = 1, 2
107 !!$ !===level_set
108 !!$ DO n = 1, inputs%nb_fluid -1
109 !!$ level_set_m1(n,:,j,i) = level_set_exact(n,j,vv_mesh%rr,mode,time-dt)
110 !!$ level_set (n,:,j,i) = level_set_exact(n,j,vv_mesh%rr,mode,time)
111 !!$ END DO
112 !!$ END DO
113 !!$ END DO
114 !!$
115 !!$ END SUBROUTINE init_level_set
116 
117 !!$ !===Source in momemtum equation. Always called.
118 !!$ FUNCTION source_in_NS_momentum(TYPE, rr, mode, i, time, Re, ty, opt_density, opt_tempn) RESULT(vv)
119 !!$ IMPLICIT NONE
120 !!$ INTEGER , INTENT(IN) :: TYPE
121 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
122 !!$ INTEGER , INTENT(IN) :: mode, i
123 !!$ REAL(KIND=8), INTENT(IN) :: time
124 !!$ REAL(KIND=8), INTENT(IN) :: Re
125 !!$ CHARACTER(LEN=2), INTENT(IN) :: ty
126 !!$ REAL(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: opt_density
127 !!$ REAL(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: opt_tempn
128 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
129 !!$
130 !!$ vv = 0.d0
131 !!$ CALL error_petsc('source_in_NS_momentum: should not be called for this test')
132 !!$ RETURN
133 !!$ END FUNCTION source_in_NS_momentum
134 
135 !!$ !===Extra source in temperature equation. Always called.
136 !!$ FUNCTION source_in_temperature(TYPE, rr, m, t)RESULT(vv)
137 !!$ IMPLICIT NONE
138 !!$ INTEGER , INTENT(IN) :: TYPE
139 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
140 !!$ INTEGER , INTENT(IN) :: m
141 !!$ REAL(KIND=8), INTENT(IN) :: t
142 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
143 !!$
144 !!$ vv = 0.d0
145 !!$ CALL error_petsc('source_in_temperature: should not be called for this test')
146 !!$ RETURN
147 !!$ END FUNCTION source_in_temperature
148 
149 !!$ !===Extra source in level set equation. Always called.
150 !!$ FUNCTION source_in_level_set(interface_nb,TYPE, rr, m, t)RESULT(vv)
151 !!$ IMPLICIT NONE
152 !!$ INTEGER , INTENT(IN) :: TYPE
153 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
154 !!$ INTEGER , INTENT(IN) :: m, interface_nb
155 !!$ REAL(KIND=8), INTENT(IN) :: t
156 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
157 !!$
158 !!$ vv=0.d0
159 !!$ CALL error_petsc('sourece_in_temperature: should not be called for this test')
160 !!$ END FUNCTION source_in_level_set
161 
162 !!$ !===Velocity for boundary conditions in Navier-Stokes.
163 !!$ !===Can be used also to initialize velocity in: init_velocity_pressure_temperature
164 !!$ FUNCTION vv_exact(TYPE,rr,m,t) RESULT(vv)
165 !!$ IMPLICIT NONE
166 !!$ INTEGER , INTENT(IN) :: TYPE
167 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
168 !!$ INTEGER, INTENT(IN) :: m
169 !!$ REAL(KIND=8), INTENT(IN) :: t
170 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
171 !!$
172 !!$ vv(:) = 0.d0
173 !!$ CALL error_petsc('vv_exact: should not be called for this test')
174 !!$ RETURN
175 !!$ END FUNCTION vv_exact
176 
177 !!$ !===Solid velocity imposed when using penalty technique
178 !!$ !===Defined in Fourier space on mode 0 only.
179 !!$ FUNCTION imposed_velocity_by_penalty(rr,t) RESULT(vv)
180 !!$ IMPLICIT NONE
181 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
182 !!$ REAL(KIND=8), INTENT(IN) :: t
183 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2),6) :: vv
184 !!$
185 !!$ vv=0.d0
186 !!$ RETURN
187 !!$ END FUNCTION imposed_velocity_by_penalty
188 
189 !!$ !===Pressure for boundary conditions in Navier-Stokes.
190 !!$ !===Can be used also to initialize pressure in the subroutine init_velocity_pressure.
191 !!$ !===Use this routine for outflow BCs only.
192 !!$ !===CAUTION: Do not enfore BCs on pressure where normal component
193 !!$ ! of velocity is prescribed.
194 !!$ FUNCTION pp_exact(TYPE,rr,m,t) RESULT (vv)
195 !!$ IMPLICIT NONE
196 !!$ INTEGER , INTENT(IN) :: TYPE
197 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
198 !!$ INTEGER , INTENT(IN) :: m
199 !!$ REAL(KIND=8), INTENT(IN) :: t
200 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
201 !!$
202 !!$ vv=0.d0
203 !!$ CALL error_petsc('pp_exact: should not be called for this test')
204 !!$ RETURN
205 !!$ END FUNCTION pp_exact
206 
207 !!$ !===Temperature for boundary conditions in temperature equation.
208 !!$ FUNCTION temperature_exact(TYPE,rr,m,t) RESULT (vv)
209 !!$ IMPLICIT NONE
210 !!$ INTEGER , INTENT(IN) :: TYPE
211 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
212 !!$ INTEGER , INTENT(IN) :: m
213 !!$ REAL(KIND=8), INTENT(IN) :: t
214 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
215 !!$
216 !!$ vv = 0.d0
217 !!$ CALL error_petsc('temperature_exact: should not be called for this test')
218 !!$ RETURN
219 !!$ END FUNCTION temperature_exact
220 
221 !!$ !===Can be used to initialize level set in the subroutine init_level_set.
222 !!$ FUNCTION level_set_exact(interface_nb,TYPE,rr,m,t) RESULT (vv)
223 !!$ IMPLICIT NONE
224 !!$ INTEGER , INTENT(IN) :: TYPE
225 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
226 !!$ INTEGER , INTENT(IN) :: m, interface_nb
227 !!$ REAL(KIND=8), INTENT(IN) :: t
228 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
229 !!$
230 !!$ vv = 0.d0
231 !!$ CALL error_petsc('level_set_exact: should not be called for this test')
232 !!$ RETURN
233 !!$
234 !!$ END FUNCTION level_set_exact
235 
236 !!$ !===Penalty coefficient (if needed)
237 !!$ !===This coefficient is equal to zero in subdomain
238 !!$ !===where penalty is applied (penalty is zero in solid)
239 !!$ FUNCTION penal_in_real_space(mesh,rr_gauss,angles,nb_angles,nb,ne,time) RESULT(vv)
240 !!$ IMPLICIT NONE
241 !!$ TYPE(mesh_type) :: mesh
242 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr_gauss
243 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: angles
244 !!$ INTEGER, INTENT(IN) :: nb_angles
245 !!$ INTEGER, INTENT(IN) :: nb, ne
246 !!$ REAL(KIND=8), INTENT(IN) :: time
247 !!$ REAL(KIND=8), DIMENSION(nb_angles,ne-nb+1) :: vv
248 !!$
249 !!$ vv = 1.d0
250 !!$ CALL error_petsc('penal_in_real_space: should not be called for this test')
251 !!$ RETURN
252 !!$ END FUNCTION penal_in_real_space
253 
254 !!$ !===Extension of the velocity field in the solid.
255 !!$ !===Used when temperature or Maxwell equations are solved.
256 !!$ !===It extends the velocity field on the Navier-Stokes domain to a
257 !!$ !===velocity field on the temperature and the Maxwell domain.
258 !!$ !===It is also used if problem type=mxw and restart velocity
259 !!$ !===is set to true in data (type problem denoted mxx in the code).
260 !!$ FUNCTION extension_velocity(TYPE, H_mesh, mode, t, n_start) RESULT(vv)
261 !!$ IMPLICIT NONE
262 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
263 !!$ INTEGER , INTENT(IN) :: TYPE, n_start
264 !!$ INTEGER, INTENT(IN) :: mode
265 !!$ REAL(KIND=8), INTENT(IN) :: t
266 !!$ REAL(KIND=8), DIMENSION(H_Mesh%np) :: vv
267 !!$
268 !!$ vv = 0.d0
269 !!$ RETURN
270 !!$
271 !!$ END FUNCTION extension_velocity
272 
273  !===============================================================================
274  ! Boundary conditions for Maxwell
275  !===============================================================================
276  !===Velocity used in the induction equation.
277  !===Used only if problem type is mxw and restart velocity is false
278  FUNCTION vexact(m, H_mesh) RESULT(vv) !Set uniquement a l'induction
279  IMPLICIT NONE
280  TYPE(mesh_type), INTENT(IN) :: h_mesh
281  INTEGER, INTENT(IN) :: m
282  REAL(KIND=8), DIMENSION(H_mesh%np,6) :: vv
283  REAL(KIND=8), DIMENSION(:), POINTER :: r, z
284 
285  IF (m==0) THEN
286  vv = 0
287  RETURN
288  END IF
289  r => h_mesh%rr(1,:)
290  z => h_mesh%rr(2,:)
291  vv(:,1) = alpha*z*(r**(m-1))*m !-alpha*z*gamma/r**(m+1)*m
292  vv(:,2) = beta *z*(r**(m-1))*m !-beta *z*gamma/r**(m+1)*m
293  vv(:,3) = beta *z*(r**(m-1))*m !+beta *z*gamma/r**(m+1)*m
294  vv(:,4) =-alpha*z*(r**(m-1))*m !-alpha*z*gamma/r**(m+1)*m
295  vv(:,5) = alpha*(r**m) !+ alpha*gamma/r**m)
296  vv(:,6) = beta *(r**m) !+ beta*gamma/r**m)
297  vv = vv/m**3
298  RETURN
299  END FUNCTION vexact
300 
301 !!$ !===Magnetic field and magnetic induction for quasi-static approximation
302 !!$ !===if needed
303 !!$ FUNCTION H_B_quasi_static(char_h_b, rr, m) RESULT(vv)
304 !!$ IMPLICIT NONE
305 !!$ CHARACTER(LEN=1), INTENT(IN) :: char_h_b
306 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
307 !!$ INTEGER, INTENT(IN) :: m
308 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2),6) :: vv
309 !!$
310 !!$ vv = 0.d0
311 !!$ RETURN
312 !!$ END FUNCTION H_B_quasi_static
313 
314  !===Magnetic field for boundary conditions in the Maxwell equations.
315  FUNCTION hexact(H_mesh, TYPE, rr, m, mu_H_field, t) RESULT(vv)
316  IMPLICIT NONE
317  TYPE(mesh_type), INTENT(IN) :: h_mesh
318  INTEGER , INTENT(IN) :: type
319  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
320  INTEGER , INTENT(IN) :: m
321  REAL(KIND=8), INTENT(IN) :: t
322  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: mu_h_field
323  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
324  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: r,z
325  REAL(KIND=8) :: muh
326 
327  IF (maxval(mu_h_field) /= minval(mu_h_field)) THEN
328  CALL error_petsc(' BUG in condlim, mu not constant')
329  END IF
330  muh=mu_h_field(1)
331  r = rr(1,:)
332  z = rr(2,:)
333  IF (m==0) THEN
334  vv = 0
335  RETURN
336  END IF
337  IF (TYPE == 1) then
338  vv = alpha*z*(r**(m-1))*m !-alpha*z*gamma/r**(m+1)*m
339  ELSEIF (TYPE == 2) then
340  vv = beta *z*(r**(m-1))*m !-beta *z*gamma/r**(m+1)*m
341  ELSEIF (TYPE ==3) then
342  vv = beta *z*(r**(m-1))*m !+beta *z*gamma/r**(m+1)*m
343  ELSEIF (TYPE == 4) then
344  vv =-alpha*z*(r**(m-1))*m !+-alpha*z*gamma/r**(m+1)*m
345  ELSEIF (TYPE == 5) then
346  vv = alpha*(r**m) ! + alpha*(gamma/r**m)
347  ELSEIF (TYPE == 6) then
348  vv = beta *(r**m) ! + beta *(gamma/r**m)
349  ENDIF
350  vv = (vv/muh)*cos(t)/m**3
351  RETURN
352 
353  !===Dummies variables to avoid warning
354  r=h_mesh%rr(1,1)
355  !===Dummies variables to avoid warning
356  END FUNCTION hexact
357 
358  !===Scalar potential for boundary conditions in the Maxwell equations.
359  FUNCTION phiexact(TYPE, rr, m, mu_phi,t) RESULT(vv)
360  IMPLICIT NONE
361  INTEGER , INTENT(IN) :: type
362  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
363  INTEGER , INTENT(IN) :: m
364  REAL(KIND=8), INTENT(IN) :: mu_phi, t
365  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
366  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: r, z
367 
368  r = rr(1,:)
369  z = rr(2,:)
370  vv(:) = 0.d0
371  IF (m==0) RETURN
372  IF (TYPE == 1) then
373  vv = alpha*z*(r**m) ! + alpha*z*(gamma/r**m)
374  ELSEIF (TYPE == 2) then
375  vv = beta *z*(r**m) ! + beta *z*(gamma/r**m)
376  ENDIF
377  vv = (vv/mu_phi)*cos(t)/m**3
378  RETURN
379  END FUNCTION phiexact
380 
381  !===Current in Ohm's law. Curl(H) = sigma(E + uxB) + current
382  FUNCTION jexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t, mesh_id, opt_B_ext) RESULT(vv)
383  IMPLICIT NONE
384  INTEGER , INTENT(IN) :: type
385  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: rr
386  INTEGER , INTENT(IN) :: m
387  REAL(KIND=8), INTENT(IN) :: mu_phi, sigma, mu_h, t
388  INTEGER , INTENT(IN) :: mesh_id
389  REAL(KIND=8), DIMENSION(6), OPTIONAL,INTENT(IN) :: opt_b_ext
390  REAL(KIND=8) :: vv
391  REAL(KIND=8) :: x !dummy
392  INTEGER :: n !dummy
393 
394  vv = -sigma* eexact_gauss(type, rr, m, mu_phi, sigma, mu_h, t)
395  RETURN
396 
397  !===Dummies variables to avoid warning
398  n=mesh_id
399  IF (present(opt_b_ext)) x=opt_b_ext(1)
400  !===Dummies variables to avoid warning
401  END FUNCTION jexact_gauss
402 
403  !===Electric field for Neumann BC (cf. doc)
404  FUNCTION eexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t) RESULT(vv)
405  IMPLICIT NONE
406  INTEGER, INTENT(IN) :: type
407  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: rr
408  INTEGER, INTENT(IN) :: m
409  REAL(KIND=8), INTENT(IN) :: mu_phi, sigma, mu_h, t
410  REAL(KIND=8) :: vv
411  REAL(KIND=8) :: r, z
412 
413  r = rr(1)
414  z = rr(2)
415  vv = 0.d0
416 
417  IF (m == 0) RETURN
418  IF (TYPE == 1) then
419  vv = 0.d0
420  ELSEIF (TYPE == 2) then
421  vv = 0.d0
422  ELSEIF (TYPE ==3) then
423  vv = alpha*(-1.d0/(m+2)*r**(m+1)) ! + alpha*(1.d0/(m-2)*gamma/r**(m-1))
424  ELSEIF (TYPE == 4) then
425  vv = beta *(-1.d0/(m+2)*r**(m+1)) !+ beta *(1.d0/(m-2)*gamma/r**(m-1))
426  ELSEIF (TYPE == 5) then
427  vv = beta*z*(r**m) ! beta*z*(-gamma/r**m)
428  ELSEIF (TYPE == 6) then
429  vv =-alpha*z*(r**m) ! -alpha*z*(-gamma/r**m)
430  ENDIF
431  vv = -vv*sin(t)/m**3
432  RETURN
433 
434  !===Dummies variables to avoid warning
435  r=mu_phi; r=sigma; r=mu_h
436  !===Dummies variables to avoid warning
437  END FUNCTION eexact_gauss
438 
439  !===Initialization of magnetic field and scalar potential (if present)
440  SUBROUTINE init_maxwell(H_mesh, phi_mesh, time, dt, mu_H_field, mu_phi, &
441  list_mode, hn1, hn, phin1, phin)
442  IMPLICIT NONE
443  TYPE(mesh_type) :: h_mesh, phi_mesh
444  REAL(KIND=8), INTENT(OUT):: time
445  REAL(KIND=8), INTENT(IN) :: dt
446  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: mu_h_field
447  REAL(KIND=8), INTENT(IN) :: mu_phi
448  INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
449  REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: hn, hn1
450  REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: phin, phin1
451  INTEGER :: i, k
452 
453  time = -dt
454  DO k=1,6
455  DO i=1, SIZE(list_mode)
456  hn1(:,k,i) = hexact(h_mesh,k, h_mesh%rr, list_mode(i), mu_h_field, time)
457  IF (inputs%nb_dom_phi>0) THEN
458  IF (k<3) THEN
459  phin1(:,k,i) = phiexact(k, phi_mesh%rr, list_mode(i) , mu_phi, time)
460  ENDIF
461  END IF
462  ENDDO
463  ENDDO
464 
465  time = time + dt
466  DO k=1,6
467  DO i=1, SIZE(list_mode)
468  hn(:,k,i) = hexact(h_mesh,k, h_mesh%rr, list_mode(i), mu_h_field, time)
469  IF (inputs%nb_dom_phi>0) THEN
470  IF (k<3) THEN
471  phin(:,k,i) = phiexact(k, phi_mesh%rr, list_mode(i), mu_phi, time)
472  ENDIF
473  END IF
474  ENDDO
475  ENDDO
476  END SUBROUTINE init_maxwell
477 
478 !!$ !===Analytical permeability (if needed)
479 !!$ !===This function is not needed unless the flag
480 !!$ !=== ===Use FEM Interpolation for magnetic permeability (true/false)
481 !!$ !===is activated and set to .FALSE. in the data data file. Default is .TRUE.
482 !!$ FUNCTION mu_bar_in_fourier_space(H_mesh,nb,ne,pts,pts_ids) RESULT(vv)
483 !!$ IMPLICIT NONE
484 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
485 !!$ REAL(KIND=8), DIMENSION(ne-nb+1) :: vv
486 !!$ INTEGER, INTENT(IN) :: nb, ne
487 !!$ REAL(KIND=8),DIMENSION(2,ne-nb+1),OPTIONAL :: pts
488 !!$ INTEGER, DIMENSION(ne-nb+1), OPTIONAL :: pts_ids
489 !!$
490 !!$ vv = 1.d0
491 !!$ CALL error_petsc('mu_bar_in_fourier_space: should not be called for this test')
492 !!$ RETURN
493 !!$ END FUNCTION mu_bar_in_fourier_space
494 
495 !!$ !===Analytical mu_in_fourier_space (if needed)
496 !!$ !===This function is not needed unless the flag
497 !!$ !=== ===Use FEM Interpolation for magnetic permeability (true/false)
498 !!$ !===is activated and set to .FALSE. in the data data file. Default is .TRUE.
499 !!$ FUNCTION grad_mu_bar_in_fourier_space(pt,pt_id) RESULT(vv)
500 !!$ IMPLICIT NONE
501 !!$ REAL(KIND=8),DIMENSION(2), INTENT(in):: pt
502 !!$ INTEGER,DIMENSION(1), INTENT(in) :: pt_id
503 !!$ REAL(KIND=8),DIMENSION(2) :: vv
504 !!$
505 !!$ vv=0.d0
506 !!$ CALL error_petsc('grad_mu_bar_in_fourier_space: should not be called for this test')
507 !!$ RETURN
508 !!$ END FUNCTION grad_mu_bar_in_fourier_space
509 
510 !!$ !===Analytical permeability, mu in real space (if needed)
511 !!$ FUNCTION mu_in_real_space(H_mesh,angles,nb_angles,nb,ne,time) RESULT(vv)
512 !!$ IMPLICIT NONE
513 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
514 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: angles
515 !!$ INTEGER, INTENT(IN) :: nb_angles
516 !!$ INTEGER, INTENT(IN) :: nb, ne
517 !!$ REAL(KIND=8), INTENT(IN) :: time
518 !!$ REAL(KIND=8), DIMENSION(nb_angles,ne-nb+1) :: vv
519 !!$
520 !!$
521 !!$ vv = 1.d0
522 !!$ CALL error_petsc('mu_in_real_space: should not be called for this test')
523 !!$ RETURN
524 !!$ END FUNCTION mu_in_real_space
525 
526 END MODULE boundary_test_3
section doc_intro_frame_work_num_app Numerical approximation subsection doc_intro_fram_work_num_app_Fourier_FEM Fourier Finite element representation The SFEMaNS code uses a hybrid Fourier Finite element formulation The Fourier decomposition allows to approximate the problem’s solutions for each Fourier mode modulo nonlinear terms that are made explicit The variables are then approximated on a meridian section of the domain with a finite element method The numerical approximation of a function f $f f is written in the following generic t
Definition: doc_intro.h:199
subroutine, public init_maxwell(H_mesh, phi_mesh, time, dt, mu_H_field, mu_phi, list_mode, Hn1, Hn, phin1, phin)
real(kind=8) function, public jexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t, mesh_id, opt_B_ext)
real(kind=8) function, dimension(h_mesh%np, 6), public vexact(m, H_mesh)
real(kind=8) function, dimension(size(rr, 2)), public phiexact(TYPE, rr, m, mu_phi, t)
real(kind=8) function, public eexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t)
subroutine error_petsc(string)
Definition: my_util.f90:15
real(kind=8) function, dimension(size(rr, 2)), public hexact(H_mesh, TYPE, rr, m, mu_H_field, t)
section doc_intro_frame_work_num_app Numerical approximation subsection doc_intro_fram_work_num_app_Fourier_FEM Fourier Finite element representation The SFEMaNS code uses a hybrid Fourier Finite element formulation The Fourier decomposition allows to approximate the problem’s solutions for each Fourier mode modulo nonlinear terms that are made explicit The variables are then approximated on a meridian section of the domain with a finite element method The numerical approximation of a function f $f f is written in the following generic z
Definition: doc_intro.h:193