53 integer,
parameter,
public :: &
54 EXTERNAL_CURRENT_PARSER = 0, &
61 type(states_mxll_t),
intent(inout) :: st
62 class(space_t),
intent(in) :: space
63 class(mesh_t),
intent(in) :: mesh
64 real(real64),
intent(in) :: time
65 complex(real64),
optional,
intent(inout) :: rs_current_density_ext(:,:)
67 real(real64),
allocatable :: current(:,:,:)
71 safe_allocate(current(1:mesh%np, 1:space%dim, 1))
76 safe_deallocate_a(current)
84 type(states_mxll_t),
intent(inout) :: st
85 class(space_t),
intent(in) :: space
86 class(mesh_t),
intent(in) :: mesh
87 type(namespace_t),
intent(in) :: namespace
90 integer :: ip, il, nlines, ncols, idir, ierr
91 real(real64) :: j_vector(space%dim), dummy(space%dim), xx(1:space%dim), rr, omega
92 character(len=1024) :: tdf_expression, phase_expression
121 if (
parse_block(namespace,
'UserDefinedMaxwellExternalCurrent', blk) == 0)
then
126 st%external_current_number = nlines
127 safe_allocate(st%external_current_modus(1:nlines))
128 safe_allocate(st%external_current_string(1:space%dim, 1:nlines))
129 safe_allocate(st%external_current_amplitude(1:mesh%np, 1:space%dim, 1:nlines))
130 safe_allocate(st%external_current_td_function(1:nlines))
131 safe_allocate(st%external_current_omega(1:nlines))
132 safe_allocate(st%external_current_td_phase(1:nlines))
138 if ((ncols /= 4) .and. (ncols /= 5) .and. (ncols /= 6) .and. (ncols /= 7))
then
139 message(1) =
'Each line in the MaxwellExternalCurrent block must have'
140 message(2) =
'four, five, six or or seven columns.'
146 if (st%external_current_modus(il) == external_current_parser)
then
148 do idir = 1, space%dim
152 do idir = 1, space%dim
157 st%external_current_string(idir, il))
158 st%external_current_amplitude(ip, idir, il) = j_vector(idir)
162 st%external_current_omega(il) = omega
164 call tdf_read(st%external_current_td_function(il), namespace, trim(tdf_expression), ierr)
167 call tdf_read(st%external_current_td_phase(il), namespace, trim(phase_expression), ierr)
169 write(
message(1),
'(3A)')
'Error in the "', trim(tdf_expression),
'" field defined in the TDExternalFields block:'
170 write(
message(2),
'(3A)')
'Time-dependent phase function "', trim(phase_expression),
'" not found.'
174 call tdf_init(st%external_current_td_phase(il))
189 class(
space_t),
intent(in) :: space
190 class(
mesh_t),
intent(in) :: mesh
191 real(real64),
intent(in) :: time
192 real(real64),
intent(inout) :: current(:,:)
194 integer :: ip, jn, idir
195 real(real64) :: xx(space%dim), rr, tt, j_vector(space%dim), dummy(space%dim), amp(space%dim)
196 complex(real64) :: exp_arg
197 real(real64) :: tmp_amp, phase
204 do jn = 1, st%external_current_number
205 if (st%external_current_modus(jn) == external_current_parser)
then
207 call mesh_r(mesh, ip, rr, coords = xx)
208 do idir = 1, space%dim
211 trim(st%external_current_string(idir,jn)))
213 current(ip, :) = current(ip, :) + j_vector
217 exp_arg = st%external_current_omega(jn) * time +
tdf(st%external_current_td_phase(jn),time)
218 phase = real(
exp(-
m_zi*exp_arg), real64)
219 tmp_amp =
tdf(st%external_current_td_function(jn), time)
221 amp = st%external_current_amplitude(ip, :, jn) * tmp_amp
222 current(ip, :) = current(ip, :) + amp(:) * phase
double exp(double __x) __attribute__((__nothrow__
subroutine, public get_rs_density_ext(st, space, mesh, time, rs_current_density_ext)
integer, parameter, public external_current_td_function
subroutine, public external_current_calculation(st, space, mesh, time, current)
subroutine, public external_current_init(st, space, namespace, mesh)
real(real64), parameter, public m_zero
complex(real64), parameter, public m_zi
This module is intended to contain "only mathematical" functions and procedures.
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
pure subroutine, public mesh_r(mesh, ip, rr, origin, coords)
return the distance to the origin for a given grid point
subroutine, public messages_warning(no_lines, all_nodes, namespace)
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public parse_block_string(blk, l, c, res, convert_to_c)
integer function, public parse_block(namespace, name, blk, check_varinfo_)
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
subroutine, public build_rs_current_state(current_state, mesh, rs_current_state, ep_field, np)
subroutine, public tdf_init(f)
subroutine, public tdf_read(f, namespace, function_name, ierr)
This function initializes "f" from the TDFunctions block.
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
This module defines the unit system, used for input and output.
type(unit_system_t), public units_inp
the units systems for reading and writing
type(unit_t), public unit_one
some special units required for particular quantities
Describes mesh distribution to nodes.