KADATH SPECTRAL SOLVER

Tutorial 3 : Using Scalar Fields

The tutorial is devoted to the use of scalar fields.

Construction and assignment

Scalar fields are constructed from the Space. They can be given explicit expressions using the methods available for the various Domains (essentially radius and Cartesian-like coordinates).

Beware that Kadath is not really intended to use analytic expressions. They probably need to be relatively simple and used only as initial guesses for the numerical resolution of the system of equations.

// Assume a spherical space has been previously constructed
// Recover the number of domains if need be
//int ndom = space.get_nbr_domains() ;

Scalar func (space) ;

// Puts 1-r^2 in domain 0
func.set_domain(0) = 1 - pow(space.get_domain(0)->get_radius(),2) ;
// Puts 1/r in the other ones
for (int d=1 ; d<ndom ; d++)
       func.set_domain(d)  = 1/space.get_domain(d)->get_radius() ;

// Prints the values at each collocation point
cout << "Scalar function func " << endl ;
cout << func << endl ;

Dealing with infinite quantities

Kadath can perform standard operations using infinite quantities. For instance using 1/r up to infinity is fine, because the resulting value is zero.

However more complicated cases can lead to nan (not a number) quantities, when the library is not able to properly compute the result. In that case, one needs to pass the right value by hand, using functions like set_val_inf.

// Use the same scalar as before
// Puts r exp(-r^2/10) in the outer domains
for (int d=1 ; d<ndom ; d++)
       func.set_domain(d)  = *space.get_domain(d)->get_radius()
* exp(- space.get_domain(d)->get_radius()*space.get_domain(d)->get_radius()/10.) ;

// Nans appear at infinty
cout << "Nans at infinity" << endl ;
cout << func(ndom-1) << endl ;

// Sets the right value at infinity by hand
func.set_val_inf(0.) ;

// No more Nans
cout << "Nans are removed"  << endl ;
cout << func(ndom-1) << endl ;

Setting the basis

At this point, the scalar field contains no information about the spectral basis, thus limiting the possibilities. Another tutorial is devoted to the basis. Here one will assume that the basis is the standard one and just use :

// Sets the standard basis for scalar fields
func.std_base() ;

Accessing the values of the field

The cout function already used, prints the value of the field at the collocation points and/or its spectral coefficients, depending on which quantity has been computed.

However, it is often useful to compute the value at a given point (i.e. not a collocation one). This is done by using the function val_point, which takes a Point as an argument. The class Point just contains the absolute (i.e. typically the X,Y,Z) coordinates of the point.

// Define a Point in a three dimensional space
Point M (3) ;

// Sets its coordinates
M.set(1) = 2.93 ; // X coordinate
M.set(2) = -3.19 ; // Y coordinate
M.set(3) = 1.76 ; // Z coordinate

// Prints the values of the scalar field at M
cout << "The point M is " << M << endl ;
cout << "Field value at M " << func.val_point(M) << endl ;

Operations on the scalar field

Using the spectral coefficients, Kadath can perform various operations on the scalar field, like taking its radial derivative or multiplying it by the radius. In this latter case, the operation is performed in the coefficient space which prevents the appearance of nan quantities and gives the right value. We can force the code to compute the values at the collocation points by invoking coef_i.

// Compute the radial derivative
Scalar derr (func.der_r()) ;
cout << "Value of the radial derivative at M " << derr.val_point(M) << endl ;

// Use the multiplication by r :
Scalar multr (func.mult_r()) ;
//Force the code to compute the values at the collocation points
multr.coef_i() ;

// No nans appear at infinity
cout << "Value in the last domain" << endl ;
cout << multr(ndom-1) << endl ;

In the last example you can notice that the value at the last collocation is not exactly zero, at it should. However, this value goes to zero exponentially when the number of radial points increases (you can easily check that). It shows that, from a numerical point of view, the computation is correct.