SFEMaNS  version 4.1 (work in progress)
Reference documentation for SFEMaNS
 All Classes Files Functions Variables Groups Pages
sub_maxwell.f90
Go to the documentation of this file.
2  PUBLIC:: maxwell_decouple
3  PRIVATE
4 
5 CONTAINS
6 
7  SUBROUTINE maxwell_decouple(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, &
8  interface_h_mu, hn, bn, phin, hn1, bn1, phin1, vel, stab_in, sigma_in, &
9  r_fourier, index_fourier, mu_h_field, mu_phi, time, dt, rem, list_mode, &
10  h_phi_per, la_h, la_pmag, la_phi, la_mhd, sigma_ns, jj_v_to_h)
11 
12  USE def_type_mesh
13  USE periodic
14  USE input_data
17  IMPLICIT NONE
18  TYPE(mesh_type), INTENT(IN) :: h_mesh, phi_mesh, pmag_mesh
19  TYPE(interface_type), INTENT(IN) :: interface_h_phi, interface_h_mu
20  INTEGER, DIMENSION(:), INTENT(IN) :: list_mode
21  REAL(KIND=8), DIMENSION(:,:,:), INTENT(INOUT) :: vel
22  REAL(KIND=8), DIMENSION(H_mesh%np,6,SIZE(list_mode)), INTENT(INOUT) :: hn, hn1
23  REAL(KIND=8), DIMENSION(H_mesh%np,6,SIZE(list_mode)), INTENT(INOUT) :: bn, bn1
24  REAL(KIND=8), DIMENSION(:,:,:), INTENT(INOUT) :: phin, phin1
25  REAL(KIND=8), DIMENSION(3), INTENT(IN) :: stab_in
26  REAL(KIND=8), INTENT(IN) :: r_fourier
27  INTEGER, INTENT(IN) :: index_fourier
28  REAL(KIND=8), INTENT(IN) :: mu_phi, time, dt, rem
29  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: sigma_in, mu_h_field
30 ! TEST
31  REAL(KIND=8), DIMENSION(:,:,:), INTENT(IN) :: sigma_ns
32  INTEGER, DIMENSION(:), INTENT(IN) :: jj_v_to_h
33 ! TEST
34  !jan 29/JLG+FL/Forget about it/We replace it by H_p_phi_per/Feb 2 2010
35  TYPE(periodic_type), INTENT(IN) :: h_phi_per
36  !jan 29/JLG+FL/Forget about it/We replace it by H_p_phi_per/Feb 2 2010
37  TYPE(petsc_csr_la) :: la_h, la_pmag, la_phi, la_mhd
38 #include "petsc/finclude/petsc.h"
39  mpi_comm, DIMENSION(:), POINTER :: comm_one_d
40 
41  IF (inputs%if_maxwell_with_H) THEN
42  CALL maxwell_decouple_with_h(comm_one_d, h_mesh, pmag_mesh, phi_mesh, interface_h_phi, &
43  interface_h_mu, hn, bn, phin, hn1, bn1, phin1, vel, stab_in, sigma_in, &
44  r_fourier, index_fourier, mu_h_field, mu_phi, time, dt, rem, list_mode, &
45  h_phi_per, la_h, la_pmag, la_phi, la_mhd, sigma_ns, jj_v_to_h)
46 !!$ CALL maxwell_decouple_with_H(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, &
47 !!$ interface_H_mu, Hn, phin, Hn1, phin1, vel, stab_in, sigma_in, &
48 !!$ R_fourier, index_fourier, mu_H_field, mu_phi, time, dt, Rem, list_mode, &
49 !!$ H_phi_per, LA_H, LA_pmag, LA_phi, LA_mhd)
50  ELSE
51 ! TEST
52  CALL maxwell_decouple_with_b(comm_one_d, h_mesh, pmag_mesh, phi_mesh, interface_h_phi, &
53  interface_h_mu, hn, bn, phin, hn1, bn1, phin1, vel, stab_in, sigma_in, &
54  r_fourier, index_fourier, mu_h_field, mu_phi, time, dt, rem, list_mode, &
55  h_phi_per, la_h, la_pmag, la_phi, la_mhd, sigma_ns, jj_v_to_h)
56 !!$ CALL maxwell_decouple_with_B(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, &
57 !!$ interface_H_mu, Hn, phin, Hn1, phin1, vel, stab_in, sigma_in, &
58 !!$ R_fourier, index_fourier, mu_H_field, mu_phi, time, dt, Rem, list_mode, &
59 !!$ H_phi_per, LA_H, LA_pmag, LA_phi, LA_mhd)
60 ! TEST
61  END IF
62 
63  END SUBROUTINE maxwell_decouple
64 
65 END MODULE update_maxwell
subroutine, public maxwell_decouple_with_b(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, interface_H_mu, Hn, Bn, phin, Hn1, Bn1, phin1, vel, stab_in, sigma_in, R_fourier, index_fourier, mu_H_field, mu_phi, time, dt_in, Rem, list_mode, H_phi_per, LA_H, LA_pmag, LA_phi, LA_mhd, sigma_ns_in, jj_v_to_H)
subroutine, public maxwell_decouple_with_h(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, interface_H_mu, Hn, Bn, phin, Hn1, Bn1, phin1, vel, stab_in, sigma_in, R_fourier, index_fourier, mu_H_field, mu_phi, time, dt, Rem, list_mode, H_phi_per, LA_H, LA_pmag, LA_phi, LA_mhd, sigma_ns_in, jj_v_to_H)
subroutine, public maxwell_decouple(comm_one_d, H_mesh, pmag_mesh, phi_mesh, interface_H_phi, interface_H_mu, Hn, Bn, phin, Hn1, Bn1, phin1, vel, stab_in, sigma_in, R_fourier, index_fourier, mu_H_field, mu_phi, time, dt, Rem, list_mode, H_phi_per, LA_H, LA_pmag, LA_phi, LA_mhd, sigma_ns, jj_v_to_H)
Definition: sub_maxwell.f90:7