SFEMaNS  version 5.3
Reference documentation for SFEMaNS
condlim_test_1.f90
Go to the documentation of this file.
1 MODULE boundary_test_1
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 CONTAINS
34  !===============================================================================
35  ! Boundary conditions for Navier-Stokes
36  !===============================================================================
37 
38  !===Initialize velocity, pressure
39  SUBROUTINE init_velocity_pressure(mesh_f, mesh_c, time, dt, list_mode, &
40  un_m1, un, pn_m1, pn, phin_m1, phin)
41  IMPLICIT NONE
42  TYPE(mesh_type) :: mesh_f, mesh_c
43  REAL(KIND=8), INTENT(OUT):: time
44  REAL(KIND=8), INTENT(IN) :: dt
45  INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
46  REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: un_m1, un
47  REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: pn_m1, pn, phin_m1, phin
48  INTEGER :: mode, i, j
49  REAL(KIND=8), DIMENSION(mesh_c%np) :: pn_m2
50 
51  time = 0.d0
52  DO i= 1, SIZE(list_mode)
53  mode = list_mode(i)
54  DO j = 1, 6
55  !===velocity
56  un_m1(:,j,i) = vv_exact(j,mesh_f%rr,mode,time-dt)
57  un(:,j,i) = vv_exact(j,mesh_f%rr,mode,time)
58  END DO
59  DO j = 1, 2
60  !===pressure
61  pn_m2(:) = pp_exact(j,mesh_c%rr,mode,time-2*dt)
62  pn_m1(:,j,i) = pp_exact(j,mesh_c%rr,mode,time-dt)
63  pn(:,j,i) = pp_exact(j,mesh_c%rr,mode,time)
64  phin_m1(:,j,i) = pn_m1(:,j,i) - pn_m2(:)
65  phin(:,j,i) = pn(:,j,i) - pn_m1(:,j,i)
66  ENDDO
67  ENDDO
68  END SUBROUTINE init_velocity_pressure
69 
70 !!$ !===Initialize temperature
71 !!$ SUBROUTINE init_temperature(mesh, time, dt, list_mode, tempn_m1, tempn)
72 !!$ IMPLICIT NONE
73 !!$ TYPE(mesh_type) :: mesh
74 !!$ REAL(KIND=8), INTENT(OUT):: time
75 !!$ REAL(KIND=8), INTENT(IN) :: dt
76 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
77 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: tempn_m1, tempn
78 !!$ INTEGER :: mode, i, j
79 !!$
80 !!$ time = 0.d0
81 !!$ DO i= 1, SIZE(list_mode)
82 !!$ mode = list_mode(i)
83 !!$ DO j = 1, 2
84 !!$ tempn_m1(:,j,i) = temperature_exact(j, mesh%rr, mode, time-dt)
85 !!$ tempn (:,j,i) = temperature_exact(j, mesh%rr, mode, time)
86 !!$ ENDDO
87 !!$ ENDDO
88 !!$ END SUBROUTINE init_temperature
89 
90 !!$ !===Initialize level_set
91 !!$ SUBROUTINE init_level_set(vv_mesh, time, &
92 !!$ dt, list_mode, level_set_m1, level_set)
93 !!$ IMPLICIT NONE
94 !!$ TYPE(mesh_type) :: vv_mesh
95 !!$ REAL(KIND=8), INTENT(OUT):: time
96 !!$ REAL(KIND=8), INTENT(IN) :: dt
97 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
98 !!$ REAL(KIND=8), DIMENSION(:,:,:,:), INTENT(OUT):: level_set, level_set_m1
99 !!$ INTEGER :: mode, i, j, n
100 !!$
101 !!$ time = 0.d0
102 !!$ DO i= 1, SIZE(list_mode)
103 !!$ mode = list_mode(i)
104 !!$ DO j = 1, 2
105 !!$ !===level_set
106 !!$ DO n = 1, inputs%nb_fluid -1
107 !!$ level_set_m1(n,:,j,i) = level_set_exact(n,j,vv_mesh%rr,mode,time-dt)
108 !!$ level_set (n,:,j,i) = level_set_exact(n,j,vv_mesh%rr,mode,time)
109 !!$ END DO
110 !!$ END DO
111 !!$ END DO
112 !!$
113 !!$ END SUBROUTINE init_level_set
114 
115  !===Source in momemtum equation. Always called.
116  FUNCTION source_in_ns_momentum(TYPE, rr, mode, i, time, Re, ty, opt_density, opt_tempn) RESULT(vv)
117  IMPLICIT NONE
118  INTEGER , INTENT(IN) :: TYPE
119  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
120  INTEGER , INTENT(IN) :: mode, i
121  REAL(KIND=8), INTENT(IN) :: time
122  REAL(KIND=8), INTENT(IN) :: Re
123  CHARACTER(LEN=2), INTENT(IN) :: ty
124  REAL(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: opt_density
125  REAL(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: opt_tempn
126  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
127  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: LAP, TEMP,GPRESS
128  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: r, z
129  REAL(KIND=8), DIMENSION(SIZE(rr,2),6) :: Rot, V, vt
130  INTEGER :: m, m1
131  REAL(KIND=8) :: t
132  CHARACTER(LEN=2) :: np
133 
134  IF (PRESENT(opt_density)) CALL error_petsc('density should not be present for test 1')
135  IF (PRESENT(opt_tempn)) CALL error_petsc('temperature should not be present for test 1')
136 
137  r = rr(1,:)
138  z = rr(2,:)
139  m=mode
140  t=time
141  IF (m==0 .OR. m==2) THEN
142  !Nonlinear contribution
143  m1 = 1
144  rot(:,1) = m1*r*z**3*(4+m1) &
145  -6*r**3*z+9*(r*z)**2
146  rot(:,2) =-m1*r*z**3*(4-m1)&
147  -6*r**3*z+9*(r*z)**2
148  rot(:,3) = 3*m1*(r*z)**2-6*r**3*z &
149  -2*r*z**3*(4-m1)
150  rot(:,4) =-3*m1*(r*z)**2-6*r**3*z &
151  -2*r*z**3*(4+m1)
152  rot(:,5) = 12*(r*z)**2-9*r*z**3 + &
153  m1*(m1*r*z**3+3*(r*z)**2)
154  rot(:,6) = 12*(r*z)**2-9*r*z**3 + &
155  m1*(m1*r*z**3-3*(r*z)**2)
156  rot = rot*cos(t)
157 
158  vt(:,1) = m1*z**3 - 3*r*z**2
159  vt(:,4) =3*(r-z)*z**2
160  vt(:,5) =(4-m1)*z**3
161  v(:,1) = vt(:,1)*r**2
162  v(:,4) = vt(:,4)*r**2
163  v(:,5) = vt(:,5)*r**2
164 
165  vt(:,2) =-m1*z**3 - 3*r*z**2
166  vt(:,3) =3*(r-z)*z**2
167  vt(:,6) =(4+m1)*z**3
168  v(:,2) = vt(:,2)*r**2
169  v(:,3) = vt(:,3)*r**2
170  v(:,6) = vt(:,6)*r**2
171  v = v*cos(t)
172 
173  IF (m==0) THEN ! vv = -VxRot(v) pour mode 0
174  IF (mod(TYPE,2)==0) then
175  vv = 0.d0
176  ELSEIF (TYPE == 1) then
177  vv = -0.5d0*(v(:,3)*rot(:,5)-v(:,5)*rot(:,3) &
178  + v(:,4)*rot(:,6)-v(:,6)*rot(:,4))
179  ELSEIF (TYPE == 3) then
180  vv = -0.5d0*(v(:,5)*rot(:,1)-v(:,1)*rot(:,5) &
181  + v(:,6)*rot(:,2)-v(:,2)*rot(:,6))
182  ELSEIF (TYPE == 5) then
183  vv = -0.5d0*(v(:,1)*rot(:,3)-v(:,3)*rot(:,1) &
184  + v(:,2)*rot(:,4)-v(:,4)*rot(:,2))
185  ENDIF
186  RETURN
187  END IF
188 
189  IF (m==2) THEN ! vv = -VxRot(v) pour mode 2
190  IF (TYPE == 1) then
191  vv = -0.5d0*(v(:,3)*rot(:,5)-v(:,5)*rot(:,3) &
192  - (v(:,4)*rot(:,6)-v(:,6)*rot(:,4)))
193  ELSE IF (TYPE == 2) then
194  vv = -0.5d0*(v(:,3)*rot(:,6)-v(:,5)*rot(:,4) &
195  + (v(:,4)*rot(:,5)-v(:,6)*rot(:,3)))
196  ELSEIF (TYPE == 3) then
197  vv = -0.5d0*(v(:,5)*rot(:,1)-v(:,1)*rot(:,5) &
198  - (v(:,6)*rot(:,2)-v(:,2)*rot(:,6)))
199  ELSEIF (TYPE == 4) then
200  vv = -0.5d0*(v(:,5)*rot(:,2)-v(:,1)*rot(:,6) &
201  + (v(:,6)*rot(:,1)-v(:,2)*rot(:,5)))
202  ELSEIF (TYPE == 5) then
203  vv = -0.5d0*(v(:,1)*rot(:,3)-v(:,3)*rot(:,1) &
204  - (v(:,2)*rot(:,4)-v(:,4)*rot(:,2)))
205  ELSEIF (TYPE == 6) then
206  vv = -0.5d0*(v(:,1)*rot(:,4)-v(:,3)*rot(:,2) &
207  + (v(:,2)*rot(:,3)-v(:,4)*rot(:,1)))
208  END IF
209 
210  RETURN
211  END IF
212 
213  END IF
214  IF (m/=1) THEN
215  vv = 0.d0
216  RETURN
217  END IF
218 
219  !Calculation for mode 1
220  vt(:,1) = m*z**3 - 3*r*z**2
221  vt(:,4) =3*(r-z)*z**2
222  vt(:,5) =(4-m)*z**3
223  v(:,1) = vt(:,1)*r**2
224  v(:,4) = vt(:,4)*r**2
225  v(:,5) = vt(:,5)*r**2
226 
227  vt(:,2) =-m*z**3 - 3*r*z**2
228  vt(:,3) =3*(r-z)*z**2
229  vt(:,6) =(4+m)*z**3
230  v(:,2) = vt(:,2)*r**2
231  v(:,3) = vt(:,3)*r**2
232  v(:,6) = vt(:,6)*r**2
233 
234  IF (TYPE == 1) then
235  lap(:) = m*4*z**3 - 27*r*z**2 + m*6*r**2*z - 6*r**3 &
236  -m**2*vt(:,1) - 2*m*vt(:,4) - vt(:,1)
237  temp(:) = v(:,1)
238  gpress(:) = 2*r*z**3*cos(t)
239  ELSEIF (TYPE == 2) then
240  lap(:) =-m*4*z**3 - 27*r*z**2 - m*6*r**2*z - 6*r**3 &
241  -m**2*vt(:,2) + 2*m*vt(:,3) - vt(:,2)
242  temp(:) = v(:,2)
243  gpress(:) = 2*r*z**3*cos(t)
244  ELSEIF (TYPE == 3) then
245  lap(:) = 3*(9*r*z**2 - 4*z**3 + 2*r**3- 6*r**2*z) &
246  -m**2*vt(:,3) + 2*m*vt(:,2) -vt(:,3)
247  temp(:) = v(:,3)
248  gpress(:) = m*r*z**3*cos(t)
249  ELSEIF (TYPE == 4) then
250  lap(:) = 3*(9*r*z**2 - 4*z**3 + 2*r**3- 6*r**2*z) &
251  -m**2*vt(:,4) -2*m*vt(:,1) -vt(:,4)
252  temp(:) = v(:,4)
253  gpress(:) = -m*r*z**3*cos(t)
254  ELSEIF (TYPE == 5) then
255  lap(:) = (4-m)*(4*z**3 + 6*r**2*z) -m**2*vt(:,5)
256  temp(:) = v(:,5)
257  gpress(:) = 3*r**2*z**2*cos(t)
258  ELSEIF (TYPE == 6) then
259  lap(:) = (4+m)*(4*z**3 + 6*r**2*z) -m**2*vt(:,6)
260  temp(:) = v(:,6)
261  gpress(:) = 3*r**2*z**2*cos(t)
262  ENDIF
263  vv(:) = -sin(t)*temp + cos(t)*(-lap/re) + gpress
264  RETURN
265 
266  !===Dummies variables to avoid warning
267  m=i; np=ty
268  !===Dummies variables to avoid warning
269  END FUNCTION source_in_ns_momentum
270 
271 !!$ !===Extra source in temperature equation. Always called.
272 !!$ FUNCTION source_in_temperature(TYPE, rr, m, t)RESULT(vv)
273 !!$ IMPLICIT NONE
274 !!$ INTEGER , INTENT(IN) :: TYPE
275 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
276 !!$ INTEGER , INTENT(IN) :: m
277 !!$ REAL(KIND=8), INTENT(IN) :: t
278 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
279 !!$
280 !!$ vv = 0.d0
281 !!$ CALL error_petsc('source_in_temperature: should not be called for this test')
282 !!$ RETURN
283 !!$ END FUNCTION source_in_temperature
284 
285 !!$ !===Extra source in level set equation. Always called.
286 !!$ FUNCTION source_in_level_set(interface_nb,TYPE, rr, m, t)RESULT(vv)
287 !!$ IMPLICIT NONE
288 !!$ INTEGER , INTENT(IN) :: TYPE
289 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
290 !!$ INTEGER , INTENT(IN) :: m, interface_nb
291 !!$ REAL(KIND=8), INTENT(IN) :: t
292 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
293 !!$
294 !!$ vv=0.d0
295 !!$ CALL error_petsc('source_in_level_set: should not be called for this test')
296 !!$
297 !!$ END FUNCTION source_in_level_set
298 
299  !===Velocity for boundary conditions in Navier-Stokes.
300  !===Can be used also to initialize velocity in: init_velocity_pressure_temperature
301  FUNCTION vv_exact(TYPE,rr,m,t) RESULT(vv)
302  IMPLICIT NONE
303  INTEGER , INTENT(IN) :: type
304  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
305  INTEGER, INTENT(IN) :: m
306  REAL(KIND=8), INTENT(IN) :: t
307  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
308  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: r, z
309 
310  r = rr(1,:)
311  z = rr(2,:)
312 
313  IF (m/=1) THEN
314  vv = 0.d0
315  RETURN
316  END IF
317 
318  IF (TYPE == 1) then
319  vv(:) = m*r**2*z**3 - 3*r**3*z**2
320  ELSEIF (TYPE == 2) then
321  vv(:) = -m*r**2*z**3 - 3*r**3*z**2
322  ELSEIF (TYPE == 3) then
323  vv(:) = 3*r**3*z**2 - 3 *r**2*z**3
324  ELSEIF (TYPE == 4) then
325  vv(:) = 3*r**3*z**2 - 3 *r**2*z**3
326  ELSEIF (TYPE == 5) then
327  vv(:) = r**2*z**3*(4-m)
328  ELSEIF (TYPE == 6) then
329  vv(:) = r**2*z**3*(4+m)
330  ENDIF
331 
332  vv(:) = vv(:) * cos(t)
333 
334  RETURN
335  END FUNCTION vv_exact
336 
337 !!$ !===Solid velocity imposed when using penalty technique
338 !!$ !===Defined in Fourier space on mode 0 only.
339 !!$ FUNCTION imposed_velocity_by_penalty(rr,t) RESULT(vv)
340 !!$ IMPLICIT NONE
341 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
342 !!$ REAL(KIND=8), INTENT(IN) :: t
343 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2),6) :: vv
344 !!$
345 !!$ vv=0
346 !!$ RETURN
347 !!$ END FUNCTION imposed_velocity_by_penalty
348 
349  !===Pressure for boundary conditions in Navier-Stokes.
350  !===Can be used also to initialize pressure in the subroutine init_velocity_pressure.
351  !===Use this routine for outflow BCs only.
352  !===CAUTION: Do not enfore BCs on pressure where normal component
353  ! of velocity is prescribed.
354  FUNCTION pp_exact(TYPE,rr,m,t) RESULT (vv)
355  IMPLICIT NONE
356  INTEGER , INTENT(IN) :: TYPE
357  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
358  INTEGER , INTENT(IN) :: m
359  REAL(KIND=8), INTENT(IN) :: t
360  REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
361 
362  IF (m/=1) THEN
363  vv = 0.d0
364  RETURN
365  END IF
366  IF (type==1.OR.type==2) THEN
367  vv(:) = rr(1,:)**2*rr(2,:)**3*cos(t)
368  ELSE
369  CALL error_petsc('pp_exact: TYPE should be equal to 1 or 2')
370  END IF
371  RETURN
372  END FUNCTION pp_exact
373 
374 !!$ !===Temperature for boundary conditions in temperature equation.
375 !!$ FUNCTION temperature_exact(TYPE,rr,m,t) RESULT (vv)
376 !!$ IMPLICIT NONE
377 !!$ INTEGER , INTENT(IN) :: TYPE
378 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
379 !!$ INTEGER , INTENT(IN) :: m
380 !!$ REAL(KIND=8), INTENT(IN) :: t
381 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
382 !!$
383 !!$ vv = 0.d0
384 !!$ CALL error_petsc('temperature_exact: should not be called for this test')
385 !!$ RETURN
386 !!$ END FUNCTION temperature_exact
387 
388 !!$ !===Can be used to initialize level set in the subroutine init_level_set.
389 !!$ FUNCTION level_set_exact(interface_nb,TYPE,rr,m,t) RESULT (vv)
390 !!$ IMPLICIT NONE
391 !!$ INTEGER , INTENT(IN) :: TYPE
392 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
393 !!$ INTEGER , INTENT(IN) :: m, interface_nb
394 !!$ REAL(KIND=8), INTENT(IN) :: t
395 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
396 !!$
397 !!$ vv = 0.d0
398 !!$ CALL error_petsc('level_set_exact: should not be called for this test')
399 !!$ RETURN
400 !!$
401 !!$ END FUNCTION level_set_exact
402 
403 !!$ !===Penalty coefficient (if needed)
404 !!$ !===This coefficient is equal to zero in subdomain
405 !!$ !===where penalty is applied (penalty is zero in solid)
406 !!$ FUNCTION penal_in_real_space(mesh,rr_gauss,angles,nb_angles,nb,ne,time) RESULT(vv)
407 !!$ IMPLICIT NONE
408 !!$ TYPE(mesh_type) :: mesh
409 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr_gauss
410 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: angles
411 !!$ INTEGER, INTENT(IN) :: nb_angles
412 !!$ INTEGER, INTENT(IN) :: nb, ne
413 !!$ REAL(KIND=8), INTENT(IN) :: time
414 !!$ REAL(KIND=8), DIMENSION(nb_angles,ne-nb+1) :: vv
415 !!$
416 !!$ vv = 1.d0
417 !!$ CALL error_petsc('penal_in_real_space: should not be called for this test')
418 !!$ RETURN
419 !!$ END FUNCTION penal_in_real_space
420 
421 !!$ !===Extension of the velocity field in the solid.
422 !!$ !===Used when temperature or Maxwell equations are solved.
423 !!$ !===It extends the velocity field on the Navier-Stokes domain to a
424 !!$ !===velocity field on the temperature and the Maxwell domain.
425 !!$ !===It is also used if problem type=mxw and restart velocity
426 !!$ !===is set to true in data (type problem denoted mxx in the code).
427 !!$ FUNCTION extension_velocity(TYPE, H_mesh, mode, t, n_start) RESULT(vv)
428 !!$ IMPLICIT NONE
429 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
430 !!$ INTEGER , INTENT(IN) :: TYPE, n_start
431 !!$ INTEGER, INTENT(IN) :: mode
432 !!$ REAL(KIND=8), INTENT(IN) :: t
433 !!$ REAL(KIND=8), DIMENSION(H_Mesh%np) :: vv
434 !!$
435 !!$ vv = 0.d0
436 !!$ RETURN
437 !!$
438 !!$ END FUNCTION extension_velocity
439 
440  !===============================================================================
441  ! Boundary conditions for Maxwell
442  !===============================================================================
443 !!$ !===Velocity used in the induction equation.
444 !!$ !===Used only if problem type is mxw and restart velocity is false
445 !!$ FUNCTION Vexact(m, H_mesh) RESULT(vv) !Set uniquement a l'induction
446 !!$ IMPLICIT NONE
447 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
448 !!$ INTEGER, INTENT(IN) :: m
449 !!$ REAL(KIND=8), DIMENSION(H_mesh%np,6) :: vv
450 !!$
451 !!$ vv = 0.d0
452 !!$ CALL error_petsc('Vexact: should not be called for this test')
453 !!$ END FUNCTION Vexact
454 
455 !!$ !===Magnetic field and magnetic induction for quasi-static approximation
456 !!$ !===if needed
457 !!$ FUNCTION H_B_quasi_static(char_h_b, rr, m) RESULT(vv)
458 !!$ IMPLICIT NONE
459 !!$ CHARACTER(LEN=1), INTENT(IN) :: char_h_b
460 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
461 !!$ INTEGER, INTENT(IN) :: m
462 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2),6) :: vv
463 !!$
464 !!$ vv = 0.d0
465 !!$ RETURN
466 !!$ END FUNCTION H_B_quasi_static
467 
468 !!$ !===Magnetic field for boundary conditions in the Maxwell equations.
469 !!$ FUNCTION Hexact(H_mesh, TYPE, rr, m, mu_H_field, t) RESULT(vv)
470 !!$ IMPLICIT NONE
471 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
472 !!$ INTEGER , INTENT(IN) :: TYPE
473 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
474 !!$ INTEGER , INTENT(IN) :: m
475 !!$ REAL(KIND=8), INTENT(IN) :: t
476 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: mu_H_field
477 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
478 !!$
479 !!$ vv=0.d0
480 !!$ CALL error_petsc('Hexact: should not be called for this test')
481 !!$ RETURN
482 !!$ END FUNCTION Hexact
483 
484 !!$ !===Scalar potential for boundary conditions in the Maxwell equations.
485 !!$ FUNCTION Phiexact(TYPE, rr, m, mu_phi,t) RESULT(vv)
486 !!$ IMPLICIT NONE
487 !!$ INTEGER , INTENT(IN) :: TYPE
488 !!$ REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: rr
489 !!$ INTEGER , INTENT(IN) :: m
490 !!$ REAL(KIND=8), INTENT(IN) :: mu_phi, t
491 !!$ REAL(KIND=8), DIMENSION(SIZE(rr,2)) :: vv
492 !!$
493 !!$ vv=0.d0
494 !!$ CALL error_petsc('Phiexact: should not be called for this test')
495 !!$ RETURN
496 !!$ END FUNCTION Phiexact
497 
498 !!$ !===Current in Ohm's law. Curl(H) = sigma(E + uxB) + current
499 !!$ FUNCTION Jexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t, mesh_id, opt_B_ext) RESULT(vv)
500 !!$ IMPLICIT NONE
501 !!$ INTEGER , INTENT(IN) :: TYPE
502 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: rr
503 !!$ INTEGER , INTENT(IN) :: m
504 !!$ REAL(KIND=8), INTENT(IN) :: mu_phi, sigma, mu_H, t
505 !!$ INTEGER , INTENT(IN) :: mesh_id
506 !!$ REAL(KIND=8), DIMENSION(6), OPTIONAL,INTENT(IN) :: opt_B_ext
507 !!$ REAL(KIND=8) :: vv
508 !!$
509 !!$ vv=0.d0
510 !!$ CALL error_petsc('Jexact_gauss: should not be called for this test')
511 !!$ RETURN
512 !!$ END FUNCTION Jexact_gauss
513 
514 !!$ !===Electric field for Neumann BC (cf. doc)
515 !!$ FUNCTION Eexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t) RESULT(vv)
516 !!$ IMPLICIT NONE
517 !!$ INTEGER, INTENT(IN) :: TYPE
518 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: rr
519 !!$ INTEGER, INTENT(IN) :: m
520 !!$ REAL(KIND=8), INTENT(IN) :: mu_phi, sigma, mu_H, t
521 !!$ REAL(KIND=8) :: vv
522 !!$
523 !!$ vv = 0.d0
524 !!$ CALL error_petsc('Eexact_gauss: should not be called for this test')
525 !!$ END FUNCTION Eexact_gauss
526 
527 !!$ !===Initialization of magnetic field and scalar potential (if present)
528 !!$ SUBROUTINE init_maxwell(H_mesh, phi_mesh, time, dt, mu_H_field, mu_phi, &
529 !!$ list_mode, Hn1, Hn, phin1, phin)
530 !!$ IMPLICIT NONE
531 !!$ TYPE(mesh_type) :: H_mesh, phi_mesh
532 !!$ REAL(KIND=8), INTENT(OUT):: time
533 !!$ REAL(KIND=8), INTENT(IN) :: dt
534 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: mu_H_field
535 !!$ REAL(KIND=8), INTENT(IN) :: mu_phi
536 !!$ INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
537 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: Hn, Hn1
538 !!$ REAL(KIND=8), DIMENSION(:,:,:), INTENT(OUT):: phin, phin1
539 !!$ INTEGER :: i, k
540 !!$
541 !!$ time = -dt
542 !!$ DO k=1,6
543 !!$ DO i=1, SIZE(list_mode)
544 !!$ Hn1(:,k,i) = Hexact(H_mesh,k, H_mesh%rr, list_mode(i), mu_H_field, time)
545 !!$ IF (inputs%nb_dom_phi>0) THEN
546 !!$ IF (k<3) THEN
547 !!$ phin1(:,k,i) = Phiexact(k, phi_mesh%rr, list_mode(i) , mu_phi, time)
548 !!$ ENDIF
549 !!$ END IF
550 !!$ ENDDO
551 !!$ ENDDO
552 !!$
553 !!$ time = time + dt
554 !!$ DO k=1,6
555 !!$ DO i=1, SIZE(list_mode)
556 !!$ Hn(:,k,i) = Hexact(H_mesh,k, H_mesh%rr, list_mode(i), mu_H_field, time)
557 !!$ IF (inputs%nb_dom_phi>0) THEN
558 !!$ IF (k<3) THEN
559 !!$ phin(:,k,i) = Phiexact(k, phi_mesh%rr, list_mode(i), mu_phi, time)
560 !!$ ENDIF
561 !!$ END IF
562 !!$ ENDDO
563 !!$ ENDDO
564 !!$ END SUBROUTINE init_maxwell
565 
566 !!$ !===Analytical permeability (if needed)
567 !!$ !===This function is not needed unless the flag
568 !!$ !=== ===Use FEM Interpolation for magnetic permeability (true/false)
569 !!$ !===is activated and set to .FALSE. in the data data file. Default is .TRUE.
570 !!$ FUNCTION mu_bar_in_fourier_space(H_mesh,nb,ne,pts,pts_ids) RESULT(vv)
571 !!$ IMPLICIT NONE
572 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
573 !!$ REAL(KIND=8), DIMENSION(ne-nb+1) :: vv
574 !!$ INTEGER, INTENT(IN) :: nb, ne
575 !!$ REAL(KIND=8),DIMENSION(2,ne-nb+1),OPTIONAL :: pts
576 !!$ INTEGER, DIMENSION(ne-nb+1), OPTIONAL :: pts_ids
577 !!$
578 !!$ vv = 1.d0
579 !!$ CALL error_petsc('mu_bar_in_fourier_space: should not be called for this test')
580 !!$ RETURN
581 !!$ END FUNCTION mu_bar_in_fourier_space
582 
583 !!$ !===Analytical mu_in_fourier_space (if needed)
584 !!$ !===This function is not needed unless the flag
585 !!$ !=== ===Use FEM Interpolation for magnetic permeability (true/false)
586 !!$ !===is activated and set to .FALSE. in the data data file. Default is .TRUE.
587 !!$ FUNCTION grad_mu_bar_in_fourier_space(pt,pt_id) RESULT(vv)
588 !!$ IMPLICIT NONE
589 !!$ REAL(KIND=8),DIMENSION(2), INTENT(in):: pt
590 !!$ INTEGER,DIMENSION(1), INTENT(in) :: pt_id
591 !!$ REAL(KIND=8),DIMENSION(2) :: vv
592 !!$
593 !!$ vv=0.d0
594 !!$ CALL error_petsc('grad_mu_bar_in_fourier_space: should not be called for this test')
595 !!$ RETURN
596 !!$ END FUNCTION grad_mu_bar_in_fourier_space
597 
598 !!$ !===Analytical permeability, mu in real space (if needed)
599 !!$ FUNCTION mu_in_real_space(H_mesh,angles,nb_angles,nb,ne,time) RESULT(vv)
600 !!$ IMPLICIT NONE
601 !!$ TYPE(mesh_type), INTENT(IN) :: H_mesh
602 !!$ REAL(KIND=8), DIMENSION(:), INTENT(IN) :: angles
603 !!$ INTEGER, INTENT(IN) :: nb_angles
604 !!$ INTEGER, INTENT(IN) :: nb, ne
605 !!$ REAL(KIND=8), INTENT(IN) :: time
606 !!$ REAL(KIND=8), DIMENSION(nb_angles,ne-nb+1) :: vv
607 !!$
608 !!$ vv = 1.d0
609 !!$ CALL error_petsc('mu_in_real_space: should not be called for this test')
610 !!$ RETURN
611 !!$ END FUNCTION mu_in_real_space
612 
613 END MODULE boundary_test_1
boundary_test_1::pp_exact
real(kind=8) function, dimension(size(rr, 2)), public pp_exact(TYPE, rr, m, t)
Definition: condlim_test_1.f90:355
boundary_test_1::init_velocity_pressure
subroutine, public init_velocity_pressure(mesh_f, mesh_c, time, dt, list_mode, un_m1, un, pn_m1, pn, phin_m1, phin)
Definition: condlim_test_1.f90:41
def_type_mesh::mesh_type
Definition: def_type_mesh.f90:74
my_util
Definition: my_util.f90:1
my_util::error_petsc
subroutine error_petsc(string)
Definition: my_util.f90:16
boundary_test_1
Definition: condlim_test_1.f90:1
boundary_test_1::source_in_ns_momentum
real(kind=8) function, dimension(size(rr, 2)), public source_in_ns_momentum(TYPE, rr, mode, i, time, Re, ty, opt_density, opt_tempn)
Definition: condlim_test_1.f90:117
input_data
Definition: read_my_data.f90:233
def_type_mesh
Definition: def_type_mesh.f90:4
boundary_test_1::vv_exact
real(kind=8) function, dimension(size(rr, 2)), public vv_exact(TYPE, rr, m, t)
Definition: condlim_test_1.f90:302