Discontinuous Galerkin Library
#include "dg/algorithm.h"
Loading...
Searching...
No Matches
dg::aRealTopology< real_type, Nd > Struct Template Referenceabstract

An abstract base class for Nd-dimensional dG grids. More...

Inheritance diagram for dg::aRealTopology< real_type, Nd >:
[legend]

Public Types

using value_type = real_type
 value type of abscissas and weights
 
using host_vector = thrust::host_vector<real_type>
 
using host_grid = RealGrid<real_type, Nd>
 Associated realisation.
 

Public Member Functions

unsigned shape (unsigned u=0) const
 \( n_u N_u\) the total number of points of an axis
 
host_vector abscissas (unsigned u=0) const
 Construct grid abscissas of the u axis.
 
host_vector weights (unsigned u=0) const
 Get the weights of the u axis.
 
std::array< unsigned, Nd > get_shape () const
 \( n_u N_u\) the total number of points of an axis
 
std::array< host_vector, Nd > get_abscissas () const
 Construct abscissas for all axes.
 
std::array< host_vector, Nd > get_weights () const
 Construct weights for all axes.
 
std::array< real_type, Nd > get_p () const
 Get left boundary point \( \vec p\).
 
std::array< real_type, Nd > get_q () const
 Get right boundary point \( \vec q\).
 
std::array< real_type, Nd > get_l () const
 Get grid length \( l_u = q_u - p_u\) for all axes.
 
std::array< real_type, Nd > get_h () const
 Get grid constant \( h_u = \frac{q_u - p_u}{N_u}\) for all axes.
 
std::array< unsigned, Nd > get_N () const
 Get number of cells \( N_u\) for all axes.
 
std::array< unsigned, Nd > get_n () const
 Get number of polynomial coefficients \( n_u\) for all axes.
 
std::array< dg::bc, Nd > get_bc () const
 Get boundary condition \( b_u\) for all axes.
 
real_type p (unsigned u=0) const
 Get left boundary point \( p_u\) for axis u.
 
real_type q (unsigned u=0) const
 Get right boundary point \( q_u\) for axis u.
 
real_type h (unsigned u=0) const
 Get grid constant \( h_u = \frac{q_u - p_u}{N_u}\) for axis u.
 
real_type l (unsigned u=0) const
 Get grid length \( l_u = q_u - p_u\) for axis u.
 
unsigned n (unsigned u=0) const
 Get number of polynomial coefficients \( n_u\) for axis u.
 
unsigned N (unsigned u=0) const
 Get number of cells \( N_u\) for axis u.
 
dg::bc bc (unsigned u=0) const
 Get boundary condition \( b_u\) for axis u.
 
RealGrid< real_type, 1 > grid (unsigned u) const
 Get axis u as a 1d grid.
 
RealGrid< real_type, 1 > axis (unsigned u) const
 An alias for "grid".
 
template<size_t Md = Nd>
real_type x0 () const
 Equivalent to p(0)
 
template<size_t Md = Nd>
real_type x1 () const
 Equivalent to p(1)
 
template<size_t Md = Nd>
real_type y0 () const
 Equivalent to p(2)
 
template<size_t Md = Nd>
real_type y1 () const
 Equivalent to q(0)
 
template<size_t Md = Nd>
real_type z0 () const
 Equivalent to q(1)
 
template<size_t Md = Nd>
real_type z1 () const
 Equivalent to q(2)
 
template<size_t Md = Nd>
real_type lx () const
 Equivalent to l(0)
 
template<size_t Md = Nd>
real_type ly () const
 Equivalent to l(1)
 
template<size_t Md = Nd>
real_type lz () const
 Equivalent to l(2)
 
template<size_t Md = Nd>
real_type hx () const
 Equivalent to h(0)
 
template<size_t Md = Nd>
real_type hy () const
 Equivalent to h(1)
 
template<size_t Md = Nd>
real_type hz () const
 Equivalent to h(2)
 
template<size_t Md = Nd>
unsigned nx () const
 Equivalent to n(0)
 
template<size_t Md = Nd>
unsigned ny () const
 Equivalent to n(1)
 
template<size_t Md = Nd>
unsigned nz () const
 Equivalent to n(2)
 
template<size_t Md = Nd>
unsigned Nx () const
 Equivalent to N(0)
 
template<size_t Md = Nd>
unsigned Ny () const
 Equivalent to N(1)
 
template<size_t Md = Nd>
unsigned Nz () const
 Equivalent to N(2)
 
template<size_t Md = Nd>
dg::bc bcx () const
 Equivalent to bc(0)
 
template<size_t Md = Nd>
dg::bc bcy () const
 Equivalent to bc(1)
 
template<size_t Md = Nd>
dg::bc bcz () const
 Equivalent to bc(2)
 
template<size_t Md = Nd>
RealGrid< real_type, 1 > gx () const
 Equivalent to grid(0)
 
template<size_t Md = Nd>
RealGrid< real_type, 1 > gy () const
 Equivalent to grid(1)
 
template<size_t Md = Nd>
RealGrid< real_type, 1 > gz () const
 Equivalent to grid(2)
 
std::array< unsigned, Nd > start () const
 Start coordinate in C-order for dg::file::NcHyperslab.
 
std::array< unsigned, Nd > count () const
 Count vector in C-order for dg::file::NcHyperslab.
 
template<size_t Md = Nd>
std::enable_if_t<(Md >=2), void > multiplyCellNumbers (real_type fx, real_type fy)
 Multiply the number of cells in the first two dimensions with a given factor.
 
template<size_t Md = Nd>
std::enable_if_t<(Md==1), void > set (unsigned new_n, unsigned new_Nx)
 Set n and N in a 1-dimensional grid.
 
template<size_t Md = Nd>
std::enable_if_t<(Md==2), void > set (unsigned new_n, unsigned new_Nx, unsigned new_Ny)
 Set n and N in a 2-dimensional grid.
 
template<size_t Md = Nd>
std::enable_if_t<(Md==3), void > set (unsigned new_n, unsigned new_Nx, unsigned new_Ny, unsigned new_Nz)
 Set n and N in a 3-dimensional grid.
 
void set (unsigned new_n, std::array< unsigned, Nd > new_N)
 Same as set( {new_n, new_n,...}, new_N);
 
void set_axis (unsigned coord, unsigned new_n, unsigned new_N)
 Set n and N for axis coord.
 
void set (std::array< unsigned, Nd > new_n, std::array< unsigned, Nd > new_N)
 Set the number of polynomials and cells.
 
void set_pq (std::array< real_type, Nd > new_p, std::array< real_type, Nd > new_q)
 Reset the boundaries of the grid.
 
void set_bcs (std::array< dg::bc, Nd > new_bcs)
 Reset the boundary conditions of the grid.
 
void set (std::array< real_type, Nd > new_p, std::array< real_type, Nd > new_q, std::array< unsigned, Nd > new_n, std::array< unsigned, Nd > new_N, std::array< dg::bc, Nd > new_bcs)
 Reset the entire grid.
 
unsigned size () const
 The total number of points.
 
void display (std::ostream &os=std::cout) const
 Display.
 
template<size_t Md = Nd>
std::enable_if_t<(Md==1), bool > contains (real_type x) const
 Check if the grid contains a point.
 
template<class Vector >
bool contains (const Vector &x) const
 Check if the grid contains a point.
 

Static Public Member Functions

static constexpr unsigned ndim ()
 Dimensionality == Nd.
 

Protected Member Functions

 ~aRealTopology ()=default
 disallow deletion through base class pointer
 
 aRealTopology ()=default
 default constructor
 
 aRealTopology (std::array< real_type, Nd > p, std::array< real_type, Nd > q, std::array< unsigned, Nd > n, std::array< unsigned, Nd > N, std::array< dg::bc, Nd > bcs)
 Construct a topology directly from points and dimensions.
 
 aRealTopology (const std::array< RealGrid< real_type, 1 >, Nd > &axes)
 Construct a topology as the product of 1d axes grids.
 
 aRealTopology (const aRealTopology &src)=default
 
aRealTopologyoperator= (const aRealTopology &src)=default
 
virtual void do_set (std::array< unsigned, Nd > new_n, std::array< unsigned, Nd > new_N)=0
 Set the number of polynomials and cells.
 
virtual void do_set_pq (std::array< real_type, Nd > new_p, std::array< real_type, Nd > new_q)=0
 Reset the boundaries of the grid.
 
virtual void do_set (std::array< dg::bc, Nd > new_bcs)=0
 Reset the boundary conditions of the grid.
 

Detailed Description

template<class real_type, size_t Nd>
struct dg::aRealTopology< real_type, Nd >

An abstract base class for Nd-dimensional dG grids.

This grid defines a discretization of the \(N_d\) dimensional hypercube given by

\[ [\vec p, \vec q] = [p_0, p_1] \times [p_1,q_1] \times ... \times [p_{N_d-1}, q_{N_d-1}]\]

Each axis \( [p_i, q_i] \) is discretized using \( N_i\) equidistant cells. Each cells is then further discretized using \( n_i\) Gauss-Legendre nodes. The Gauss-Legendre nodes are tabulated by the dg::DLT class for \( 1\leq n_i \leq 20\). Each axis further can have a boundary condition \( b_i \) that is given by dg::bc For more information dG methods see

See also
Introduction to dg methods

This class in essence provides a collection of getters and setters for the aforementioned parameters together with the abscissas and weights members that are necessary for dg::evaluate and dg::create::weights. Lastly, we provide start and count members such that the grid can be used as a dg::file::NcHyperslab in NetCDF output in dg::file

For code readability in many physical contexts the first 3 dimensions get special names \( x,\ y,\ z\) i.e. \( \vec p = (x_0, y_0, z_0)\) and \( \vec q = (x_1, y_1, z_1)\). The class provides coresponding getters and setters like e.g. x0() or Nx() for p(0) and N(0)

Template Parameters
real_typeDetermines value type of abscissas and weights
NdThe number of dimensions \( N_d\)

Member Typedef Documentation

◆ host_grid

template<class real_type , size_t Nd>
using dg::aRealTopology< real_type, Nd >::host_grid = RealGrid<real_type, Nd>

Associated realisation.

◆ host_vector

template<class real_type , size_t Nd>
using dg::aRealTopology< real_type, Nd >::host_vector = thrust::host_vector<real_type>

vector type of abscissas and weights; Is used to recognise shared topology (vs MPI) dg::is_vector_v< typename Topology::host_vector, SharedVectorTag>

◆ value_type

template<class real_type , size_t Nd>
using dg::aRealTopology< real_type, Nd >::value_type = real_type

value type of abscissas and weights

Constructor & Destructor Documentation

◆ ~aRealTopology()

template<class real_type , size_t Nd>
dg::aRealTopology< real_type, Nd >::~aRealTopology ( )
protecteddefault

disallow deletion through base class pointer

◆ aRealTopology() [1/4]

template<class real_type , size_t Nd>
dg::aRealTopology< real_type, Nd >::aRealTopology ( )
protecteddefault

default constructor

◆ aRealTopology() [2/4]

template<class real_type , size_t Nd>
dg::aRealTopology< real_type, Nd >::aRealTopology ( std::array< real_type, Nd > p,
std::array< real_type, Nd > q,
std::array< unsigned, Nd > n,
std::array< unsigned, Nd > N,
std::array< dg::bc, Nd > bcs )
inlineprotected

Construct a topology directly from points and dimensions.

Parameters
pleft boundary point
qright boundary point
nnumber of polynomial coefficients for each axis
Nnumber of cells for each axis
bcsboundary condition for each axis

◆ aRealTopology() [3/4]

template<class real_type , size_t Nd>
dg::aRealTopology< real_type, Nd >::aRealTopology ( const std::array< RealGrid< real_type, 1 >, Nd > & axes)
inlineprotected

Construct a topology as the product of 1d axes grids.

Parameters
axesOne-dimensional grids for each dimension

◆ aRealTopology() [4/4]

template<class real_type , size_t Nd>
dg::aRealTopology< real_type, Nd >::aRealTopology ( const aRealTopology< real_type, Nd > & src)
protecteddefault

explicit copy constructor (default)

Parameters
srcsource

Member Function Documentation

◆ abscissas()

template<class real_type , size_t Nd>
host_vector dg::aRealTopology< real_type, Nd >::abscissas ( unsigned u = 0) const
inline

Construct grid abscissas of the u axis.

Parameters
uAxis number u<Nd
Returns
Vector \( \vec a_u\) containing abscissas for axis u
See also
dg::evaluate dg::DLT

◆ axis()

template<class real_type , size_t Nd>
RealGrid< real_type, 1 > dg::aRealTopology< real_type, Nd >::axis ( unsigned u) const
inline

An alias for "grid".

Parameters
uAxis number u<Nd
Returns
One dimensional grid

◆ bc()

template<class real_type , size_t Nd>
dg::bc dg::aRealTopology< real_type, Nd >::bc ( unsigned u = 0) const
inline

Get boundary condition \( b_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ bcx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
dg::bc dg::aRealTopology< real_type, Nd >::bcx ( ) const
inline

Equivalent to bc(0)

◆ bcy()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
dg::bc dg::aRealTopology< real_type, Nd >::bcy ( ) const
inline

Equivalent to bc(1)

◆ bcz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
dg::bc dg::aRealTopology< real_type, Nd >::bcz ( ) const
inline

Equivalent to bc(2)

◆ contains() [1/2]

template<class real_type , size_t Nd>
template<class Vector >
bool dg::aRealTopology< real_type, Nd >::contains ( const Vector & x) const
inline

Check if the grid contains a point.

Used for example in integrate_in_domain method of dg::AdaptiveTimeloop

Note
doesn't check periodicity!!
Parameters
xpoint to check
Returns
true if p0[u]<=x[u]<=p1[u] for all u, false else
Attention
returns false if x[u] is NaN or INF

◆ contains() [2/2]

template<class real_type , size_t Nd>
template<size_t Md = Nd>
std::enable_if_t<(Md==1), bool > dg::aRealTopology< real_type, Nd >::contains ( real_type x) const
inline

Check if the grid contains a point.

Used for example in integrate_in_domain method of dg::AdaptiveTimeloop

Note
doesn't check periodicity!!
Parameters
xpoint to check
Returns
true if p0[u]<=x[u]<=p1[u] for all u, false else
Attention
returns false if x[u] is NaN or INF

◆ count()

template<class real_type , size_t Nd>
std::array< unsigned, Nd > dg::aRealTopology< real_type, Nd >::count ( ) const
inline

Count vector in C-order for dg::file::NcHyperslab.

Used to construct dg::file::NcHyperslab together with start()

Returns
reverse( get_shape())
Note
In C-order the fastest dimension is the last one while our dg::evaluate and dg::kronecker make the 0 dimension/ 1st argument the fastest varying, so we return the reverse order of get_shape()

◆ display()

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::display ( std::ostream & os = std::cout) const
inline

Display.

Parameters
osoutput stream

◆ do_set() [1/2]

template<class real_type , size_t Nd>
virtual void dg::aRealTopology< real_type, Nd >::do_set ( std::array< dg::bc, Nd > new_bcs)
protectedpure virtual

Reset the boundary conditions of the grid.

Parameters
new_bcsnew boundary condition

◆ do_set() [2/2]

template<class real_type , size_t Nd>
virtual void dg::aRealTopology< real_type, Nd >::do_set ( std::array< unsigned, Nd > new_n,
std::array< unsigned, Nd > new_N )
protectedpure virtual

Set the number of polynomials and cells.

Parameters
new_nnew number of Gaussian nodes in each dimension
new_Nnew number of cells in each dimension

◆ do_set_pq()

template<class real_type , size_t Nd>
virtual void dg::aRealTopology< real_type, Nd >::do_set_pq ( std::array< real_type, Nd > new_p,
std::array< real_type, Nd > new_q )
protectedpure virtual

Reset the boundaries of the grid.

Parameters
new_pnew left boundary
new_qnew right boundary ( > x0)

◆ get_abscissas()

template<class real_type , size_t Nd>
std::array< host_vector, Nd > dg::aRealTopology< real_type, Nd >::get_abscissas ( ) const
inline

Construct abscissas for all axes.

Returns
\( \vec a_u\) for all \( u < N_d\)

◆ get_bc()

template<class real_type , size_t Nd>
std::array< dg::bc, Nd > dg::aRealTopology< real_type, Nd >::get_bc ( ) const
inline

Get boundary condition \( b_u\) for all axes.

Returns
Boundary condition \( b_u\) for all u

◆ get_h()

template<class real_type , size_t Nd>
std::array< real_type, Nd > dg::aRealTopology< real_type, Nd >::get_h ( ) const
inline

Get grid constant \( h_u = \frac{q_u - p_u}{N_u}\) for all axes.

Returns
Grid constant \( h_u\) for all u

◆ get_l()

template<class real_type , size_t Nd>
std::array< real_type, Nd > dg::aRealTopology< real_type, Nd >::get_l ( ) const
inline

Get grid length \( l_u = q_u - p_u\) for all axes.

Returns
Grid length \( l_u\) for all u

◆ get_N()

template<class real_type , size_t Nd>
std::array< unsigned, Nd > dg::aRealTopology< real_type, Nd >::get_N ( ) const
inline

Get number of cells \( N_u\) for all axes.

Returns
Number of cells \( N_u\) for all u

◆ get_n()

template<class real_type , size_t Nd>
std::array< unsigned, Nd > dg::aRealTopology< real_type, Nd >::get_n ( ) const
inline

Get number of polynomial coefficients \( n_u\) for all axes.

Returns
Number of polynomial coefficients \( n_u\) for all u

◆ get_p()

template<class real_type , size_t Nd>
std::array< real_type, Nd > dg::aRealTopology< real_type, Nd >::get_p ( ) const
inline

Get left boundary point \( \vec p\).

Returns
Left boundary \( \vec p\)

◆ get_q()

template<class real_type , size_t Nd>
std::array< real_type, Nd > dg::aRealTopology< real_type, Nd >::get_q ( ) const
inline

Get right boundary point \( \vec q\).

Returns
right boundary \( \vec q\)

◆ get_shape()

template<class real_type , size_t Nd>
std::array< unsigned, Nd > dg::aRealTopology< real_type, Nd >::get_shape ( ) const
inline

\( n_u N_u\) the total number of points of an axis

Returns
\( n_u N_u\) for all \( u < N_d\)

◆ get_weights()

template<class real_type , size_t Nd>
std::array< host_vector, Nd > dg::aRealTopology< real_type, Nd >::get_weights ( ) const
inline

Construct weights for all axes.

Returns
\( \vec w_u\) for all \( u < N_d\)

◆ grid()

template<class real_type , size_t Nd>
RealGrid< real_type, 1 > dg::aRealTopology< real_type, Nd >::grid ( unsigned u) const
inline

Get axis u as a 1d grid.

Parameters
uAxis number u<Nd
Returns
One dimensional grid

◆ gx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
RealGrid< real_type, 1 > dg::aRealTopology< real_type, Nd >::gx ( ) const
inline

Equivalent to grid(0)

◆ gy()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
RealGrid< real_type, 1 > dg::aRealTopology< real_type, Nd >::gy ( ) const
inline

Equivalent to grid(1)

◆ gz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
RealGrid< real_type, 1 > dg::aRealTopology< real_type, Nd >::gz ( ) const
inline

Equivalent to grid(2)

◆ h()

template<class real_type , size_t Nd>
real_type dg::aRealTopology< real_type, Nd >::h ( unsigned u = 0) const
inline

Get grid constant \( h_u = \frac{q_u - p_u}{N_u}\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ hx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::hx ( ) const
inline

Equivalent to h(0)

◆ hy()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::hy ( ) const
inline

Equivalent to h(1)

◆ hz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::hz ( ) const
inline

Equivalent to h(2)

◆ l()

template<class real_type , size_t Nd>
real_type dg::aRealTopology< real_type, Nd >::l ( unsigned u = 0) const
inline

Get grid length \( l_u = q_u - p_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ lx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::lx ( ) const
inline

Equivalent to l(0)

◆ ly()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::ly ( ) const
inline

Equivalent to l(1)

◆ lz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::lz ( ) const
inline

Equivalent to l(2)

◆ multiplyCellNumbers()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
std::enable_if_t<(Md >=2), void > dg::aRealTopology< real_type, Nd >::multiplyCellNumbers ( real_type fx,
real_type fy )
inline

Multiply the number of cells in the first two dimensions with a given factor.

With this function you can resize the grid ignorantly of its current size

Parameters
fxnew global number of cells is fx*Nx()
fynew global number of cells is fy*Ny() The remaining dimensions are left unchanged

◆ N()

template<class real_type , size_t Nd>
unsigned dg::aRealTopology< real_type, Nd >::N ( unsigned u = 0) const
inline

Get number of cells \( N_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ n()

template<class real_type , size_t Nd>
unsigned dg::aRealTopology< real_type, Nd >::n ( unsigned u = 0) const
inline

Get number of polynomial coefficients \( n_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ ndim()

template<class real_type , size_t Nd>
static constexpr unsigned dg::aRealTopology< real_type, Nd >::ndim ( )
inlinestaticconstexpr

Dimensionality == Nd.

◆ Nx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::Nx ( ) const
inline

Equivalent to N(0)

◆ nx()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::nx ( ) const
inline

Equivalent to n(0)

◆ Ny()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::Ny ( ) const
inline

Equivalent to N(1)

◆ ny()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::ny ( ) const
inline

Equivalent to n(1)

◆ Nz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::Nz ( ) const
inline

Equivalent to N(2)

◆ nz()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
unsigned dg::aRealTopology< real_type, Nd >::nz ( ) const
inline

Equivalent to n(2)

◆ operator=()

template<class real_type , size_t Nd>
aRealTopology & dg::aRealTopology< real_type, Nd >::operator= ( const aRealTopology< real_type, Nd > & src)
protecteddefault

explicit assignment operator (default)

Parameters
srcsource

◆ p()

template<class real_type , size_t Nd>
real_type dg::aRealTopology< real_type, Nd >::p ( unsigned u = 0) const
inline

Get left boundary point \( p_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ q()

template<class real_type , size_t Nd>
real_type dg::aRealTopology< real_type, Nd >::q ( unsigned u = 0) const
inline

Get right boundary point \( q_u\) for axis u.

Parameters
uAxis number u<Nd
Returns
Value for axis u

◆ set() [1/6]

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set ( std::array< real_type, Nd > new_p,
std::array< real_type, Nd > new_q,
std::array< unsigned, Nd > new_n,
std::array< unsigned, Nd > new_N,
std::array< dg::bc, Nd > new_bcs )
inline

Reset the entire grid.

Parameters
new_pnew left boundary
new_qnew right boundary ( > x0)
new_nnew number of Gaussian nodes in each dimension
new_Nnew number of cells in each dimension
new_bcsnew boundary condition

◆ set() [2/6]

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set ( std::array< unsigned, Nd > new_n,
std::array< unsigned, Nd > new_N )
inline

Set the number of polynomials and cells.

Parameters
new_nnew number of Gaussian nodes in each dimension
new_Nnew number of cells in each dimension

◆ set() [3/6]

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set ( unsigned new_n,
std::array< unsigned, Nd > new_N )
inline

Same as set( {new_n, new_n,...}, new_N);

◆ set() [4/6]

template<class real_type , size_t Nd>
template<size_t Md = Nd>
std::enable_if_t<(Md==1), void > dg::aRealTopology< real_type, Nd >::set ( unsigned new_n,
unsigned new_Nx )
inline

Set n and N in a 1-dimensional grid.

Parameters
new_nnew number of Gaussian nodes
new_Nxnew number of cells

◆ set() [5/6]

template<class real_type , size_t Nd>
template<size_t Md = Nd>
std::enable_if_t<(Md==2), void > dg::aRealTopology< real_type, Nd >::set ( unsigned new_n,
unsigned new_Nx,
unsigned new_Ny )
inline

Set n and N in a 2-dimensional grid.

Parameters
new_nnew number of Gaussian nodes in x and y
new_Nxnew number of cells
new_Nynew number of cells

◆ set() [6/6]

template<class real_type , size_t Nd>
template<size_t Md = Nd>
std::enable_if_t<(Md==3), void > dg::aRealTopology< real_type, Nd >::set ( unsigned new_n,
unsigned new_Nx,
unsigned new_Ny,
unsigned new_Nz )
inline

Set n and N in a 3-dimensional grid.

Same as set({new_n,new_n,1}, {new_Nx,new_Ny,new_Nz})

Parameters
new_nnew number of Gaussian nodes in x and y
Attention
Set nz to 1
Parameters
new_Nxnew number of cells in x
new_Nynew number of cells in y
new_Nznew number of cells in z

◆ set_axis()

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set_axis ( unsigned coord,
unsigned new_n,
unsigned new_N )
inline

Set n and N for axis coord.

Parameters
coordAxis coord<Nd
new_nnew number of Gaussian nodes of axis coord
new_Nnew number of cells

◆ set_bcs()

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set_bcs ( std::array< dg::bc, Nd > new_bcs)
inline

Reset the boundary conditions of the grid.

Parameters
new_bcsnew boundary condition

◆ set_pq()

template<class real_type , size_t Nd>
void dg::aRealTopology< real_type, Nd >::set_pq ( std::array< real_type, Nd > new_p,
std::array< real_type, Nd > new_q )
inline

Reset the boundaries of the grid.

Parameters
new_pnew left boundary
new_qnew right boundary ( > x0)

◆ shape()

template<class real_type , size_t Nd>
unsigned dg::aRealTopology< real_type, Nd >::shape ( unsigned u = 0) const
inline

\( n_u N_u\) the total number of points of an axis

Parameters
uAxis number u<Nd
Returns
\( n_u N_u\)

◆ size()

template<class real_type , size_t Nd>
unsigned dg::aRealTopology< real_type, Nd >::size ( ) const
inline

The total number of points.

Returns
\( \prod_{i=0}^{N-1} n_i N_i\)

◆ start()

template<class real_type , size_t Nd>
std::array< unsigned, Nd > dg::aRealTopology< real_type, Nd >::start ( ) const
inline

Start coordinate in C-order for dg::file::NcHyperslab.

Used to construct dg::file::NcHyperslab together with count()

Returns
{0}

◆ weights()

template<class real_type , size_t Nd>
host_vector dg::aRealTopology< real_type, Nd >::weights ( unsigned u = 0) const
inline

Get the weights of the u axis.

Parameters
uAxis number u<Nd
Returns
Vector \( \vec w_u\) containing weights for axis u
See also
dg::create::weights dg::DLT

◆ x0()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::x0 ( ) const
inline

Equivalent to p(0)

◆ x1()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::x1 ( ) const
inline

Equivalent to p(1)

◆ y0()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::y0 ( ) const
inline

Equivalent to p(2)

◆ y1()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::y1 ( ) const
inline

Equivalent to q(0)

◆ z0()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::z0 ( ) const
inline

Equivalent to q(1)

◆ z1()

template<class real_type , size_t Nd>
template<size_t Md = Nd>
real_type dg::aRealTopology< real_type, Nd >::z1 ( ) const
inline

Equivalent to q(2)


The documentation for this struct was generated from the following file: