Function, e.g. the density are given by an array rho(i,spin), where i ranges over the mesh points associated with the given domain (see domain decomposition), representing
rho( ri, spin )
mesh_x_global(mesh, ip) ?
mesh_x_global(mesh, ip) returns a
FLOAT vector, corresponding to the coordinates of point
ip of the global mesh.
function mesh_x_global(mesh, ipg) result(xx) type(mesh_t), intent(in) :: mesh integer(i8), intent(in) :: ipg FLOAT :: xx(1:mesh%box%dim) FLOAT :: chi(1:mesh%box%dim) integer :: ix(1:mesh%box%dim) ! no push_sub because function is called too frequently call mesh_global_index_to_coords(mesh, ipg, ix) chi(1:mesh%box%dim) = ix(1:mesh%box%dim) * mesh%spacing(1:mesh%box%dim) xx = mesh%coord_system%to_cartesian(chi(1:mesh%box%dim)) end function mesh_x_global
Integration is implemented as a straightforward summation:
do ip = 1, mesh%np dd = dd + ff(ip) end do
Taking derivatives is done by finite differences. The points involved in a derivative are defined by the stencil (see below).
Derivatives are discussed in a separate document Derivatives.md.
Note on packed states:
Note on curvilinear meshes:
mesh::x(:,:) array always contains a regular mesh, which gets ‘distorded’ to a curvilinear mesh by additional function calls.