33 REAL(KIND=8) :: pi=acos(-1.d0)
278 FUNCTION vexact(m, H_mesh) RESULT(vv) !Set uniquement a l'induction
281 INTEGER,
INTENT(IN) :: m
282 REAL(KIND=8),
DIMENSION(H_mesh%np,6) :: vv
283 LOGICAL,
DIMENSION(H_mesh%np) :: virgin
284 INTEGER :: i, mjl, njl
285 REAL(Kind=8) :: rr, zz
286 REAL(KIND=8) :: eps=1.d-5, height=2.d0, zeta=30.d0, amp_mnd=1.d0, amp_fl=0.d0
291 DO mjl = 1, h_mesh%me
294 DO njl = 1, h_mesh%gauss%n_w
295 i = h_mesh%jj(njl,mjl)
296 IF (.NOT.virgin(i)) cycle
302 vv(i,1) = amp_mnd*(-(pi/height)*rr*(1.d0-rr)**2*(1.d0+2.d0*rr)*cos(2.d0*pi*zz/height))
303 vv(i,3) = amp_mnd*(4.d0*height/pi)*rr*(1.d0-rr)*asin(zz)
304 IF (zz .LE. (-eps))
THEN 305 vv(i,3) = vv(i,3)+ amp_fl*rr*sin(pi*rr)*(1.-2.*zeta*(zz+1.)**2)*exp(-zeta*(zz+1.)**2)
306 ELSE IF (zz .GE. eps)
THEN 307 vv(i,3) = vv(i,3)+ amp_fl*rr*sin(pi*rr)*(1.-2.*zeta*(zz-1.)**2)*exp(-zeta*(zz-1.)**2)
309 vv(i,5) = amp_mnd*(1.-rr)*(1.+rr-5.*rr**2)*sin(2.*pi*zz/height)
330 FUNCTION hexact(H_mesh,TYPE, rr, m, mu_H_field, t) RESULT(vv)
333 INTEGER ,
INTENT(IN) :: TYPE
334 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
335 INTEGER ,
INTENT(IN) :: m
336 REAL(KIND=8),
INTENT(IN) :: t
337 REAL(KIND=8),
DIMENSION(:),
INTENT(IN) :: mu_H_field
338 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv
346 n=h_mesh%np; n=type; n=
SIZE(rr,1); n=m; r=mu_h_field(1); r=t
351 FUNCTION phiexact(TYPE, rr, m, mu_phi,t) RESULT(vv)
353 INTEGER ,
INTENT(IN) :: TYPE
354 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
355 INTEGER ,
INTENT(IN) :: m
356 REAL(KIND=8),
INTENT(IN) :: mu_phi, t
357 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv
362 CALL error_petsc(
'Phiexact: should not be called for this test')
366 n=type; n=
SIZE(rr,1); n=m; r=mu_phi; r=t
371 FUNCTION jexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t, mesh_id, opt_B_ext) RESULT(vv)
373 INTEGER ,
INTENT(IN) :: TYPE
374 REAL(KIND=8),
DIMENSION(:),
INTENT(IN) :: rr
375 INTEGER ,
INTENT(IN) :: m
376 REAL(KIND=8),
INTENT(IN) :: mu_phi, sigma, mu_H, t
377 INTEGER ,
INTENT(IN) :: mesh_id
378 REAL(KIND=8),
DIMENSION(6),
OPTIONAL,
INTENT(IN) :: opt_B_ext
387 r=rr(1); r=mu_phi; r=sigma; r=mu_h; r=t; n=type; n=m; n=mesh_id
388 IF (
PRESENT(opt_b_ext)) r=opt_b_ext(1)
393 FUNCTION eexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t) RESULT(vv)
395 INTEGER,
INTENT(IN) :: TYPE
396 REAL(KIND=8),
DIMENSION(:),
INTENT(IN) :: rr
397 INTEGER,
INTENT(IN) :: m
398 REAL(KIND=8),
INTENT(IN) :: mu_phi, sigma, mu_H, t
407 r=rr(1); r=mu_phi; r=sigma; r=mu_h; r=t; n=type; n=m
412 SUBROUTINE init_maxwell(H_mesh, phi_mesh, time, dt, mu_H_field, mu_phi, &
413 list_mode, hn1, hn, phin1, phin)
416 REAL(KIND=8),
INTENT(OUT):: time
417 REAL(KIND=8),
INTENT(IN) :: dt
418 REAL(KIND=8),
DIMENSION(:),
INTENT(IN) :: mu_H_field
419 REAL(KIND=8),
INTENT(IN) :: mu_phi
420 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
421 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(OUT):: Hn, Hn1
422 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(OUT):: phin, phin1
424 REAL(KIND=8) :: Ri, Ro, A
435 IF (h_mesh%me /= 0)
THEN 436 DO i = 1,
SIZE(list_mode)
437 IF (list_mode(i) == 1)
THEN 438 hn1(:,1,i) = a*(sin(n0*h_mesh%rr(2,:))/h_mesh%rr(1,:))*(ri-h_mesh%rr(1,:))**2 &
439 * (ro-h_mesh%rr(1,:))**2
440 hn(:,1,i) = a*(sin(n0*h_mesh%rr(2,:))/h_mesh%rr(1,:))*(ri-h_mesh%rr(1,:))**2 &
441 * (ro-h_mesh%rr(1,:))**2
442 hn1(:,4,i) = -2*a*sin(n0*h_mesh%rr(2,:))*(h_mesh%rr(1,:)-ri)&
443 *(h_mesh%rr(1,:)-ro)*(2.d0*h_mesh%rr(1,:)-ro-ri)
444 hn(:,4,i) = -2*a*sin(n0*h_mesh%rr(2,:))*(h_mesh%rr(1,:)-ri)&
445 *(h_mesh%rr(1,:)-ro)*(2.d0*h_mesh%rr(1,:)-ro-ri)
452 i=phi_mesh%np; a=time; a=dt; a=mu_h_field(1); a=mu_phi
real(kind=8) function, public eexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t)
subroutine error_petsc(string)
real(kind=8) function, dimension(size(rr, 2)), public phiexact(TYPE, rr, m, mu_phi, t)
real(kind=8) function, public jexact_gauss(TYPE, rr, m, mu_phi, sigma, mu_H, t, mesh_id, opt_B_ext)
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, dimension(size(rr, 2)), public hexact(H_mesh, TYPE, rr, m, mu_H_field, t)
real(kind=8) function, dimension(h_mesh%np, 6), public vexact(m, H_mesh)