Vectors and Tensors

Tolosa-lib contains different structures for vectors and tensors in order to enhance numerical implementation and lisibility. Therefore, it enables one to implement numerical schemes easily.

Vectors

Two structures for 2D and 3D vectors are implemented. The implementation of various operations is enhanced. One can easily add and substract vectors, multiply and divide vectors by a scalar, compute the euclidean norm and its square, and the product and cross products of two vectors.

2D

The structure of a 2D vector contains two components x and y and several vector procedures. This structure overloads the assignment operator and other operators by creating a pointer component for each of them.

!===================================================================================================================!
!  2D Vector Structure
!===================================================================================================================!

TYPE vec2d

  real(rp) :: x = 0._rp
  real(rp) :: y = 0._rp

CONTAINS

  procedure, pass( self ) :: equal_vec2d
  procedure, pass( self ) :: equal_vec2d_scal
  procedure               :: add_vec2d
  procedure               :: sub_vec2d , sub_vec2d_unary
  procedure, pass( self ) :: scalar_vec2d
  procedure, pass( self ) :: vec2d_scalar
  procedure               :: vec2d_div_scalar
  procedure               :: norm_vec2d
  procedure               :: square_vec2d
  procedure               :: dot_product_vec2d
  procedure               :: tens_prod2d
  procedure               :: rot_vec2d
  procedure               :: invrot_vec2d

  generic :: assignment(=) => equal_vec2d , equal_vec2d_scal

  generic :: operator(+) => add_vec2d
  generic :: operator(-) => sub_vec2d , sub_vec2d_unary
  generic :: operator(*) => scalar_vec2d , vec2d_scalar
  generic :: operator(/) => vec2d_div_scalar

  generic :: operator(.norm.    ) => norm_vec2d
  generic :: operator(.square.  ) => square_vec2d
  generic :: operator(.dotprod. ) => dot_product_vec2d
  generic :: operator(.tensprod.) => tens_prod2d
  generic :: operator(.rot.     ) => rot_vec2d
  generic :: operator(.invrot.  ) => invrot_vec2d

END TYPE vec2d

The assignement operator = points to the equal_vec2d and equal_vec2d_scal procedures. This operation assigns a value to each component of a vector.

Assignement =

• a = b (where a and b are two vec2d objects) assigns b%x to a%x and b%y to a%y
• a = r (where a is a vec2d object and r is a scalar) assigns r to a%x and a%y

The + operator points to the add_vec2d procedure. This operation adds the components of each vector.

Operation +

c = a + b (where a, b, c are vec2dobjects) assigns a%x + b%x to c%x and a%y + b%y to c%y.

The - operator points to the sub_vec2d and sub_vec2d_unary procedures. This operation substract the components of one vector to another, or substract the components of a vector by a scalar.

Operation -

• c = a - b (where a, b, c are vec2dobjects) assigns a%x - b%x to c%x and a%y - b%y to c%y.
• c = - a (where a, c are vec2dobjects) assigns -a%x to c%x and -a%y to c%y.

The * operator points to the scalar_vec2d and vec2d_scalar procedures. This operation multiply each component of a vector by a scalar.

Operation *

c = r * a or c = a * r (where a and c are vec2d objects, and r a scalar) assign r * a%x to c%x and r * a%y to c%y.

The / operator points to the vec2d_div_scalar procedure. This operation divides each component by a scalar.

Operation /

c = a / r (where a and c are vec2d objects and r a scalar) assigns a%x / r to c%x and a%y / r to c%y.

The .norm. operator points to the norm_vec2d procedure. This operation computes the euclidean norm of a vector.

Operation .norm.

r = .norm. a (where a is a vec2d object and r a scalar) assigns sqrt( a%x**2 + a%y**2 ) to r.

The .square. operator points to the square_vec2d procedure. This operatio computes the square of the euclidean norm of a vector.

Operation .square.

r = .square. a (where a is a vec2d object and r a scalar) assigns a%x**2 + a%y**2 to r.

The .dotprod. operator points to the dot_product_vec2d procedure. This operation computes the dot product of two vectors.

Operation .dotprod.

r = a .dotprod. b (where a and b are vec2d objects and r a scalar) assigns a%x * b%x + a%y * b%y to r.

The .tensprod. operator points to the tens_prod2d procedure. This operation computes the tensor product of two vectors.

Operation .tensprod.

c = a .tensprod. b (where a, b are vec2d objects and c a tens2d object) :

• c%xx = a%x * b%x
• c%xy = a%x * b%y
• c%yx = a%y * b%x
• c%yy = a%y * b%y

The .rot. operator points to the rot_vec2d procedure. This operation applies a rotation operation on a vector.

Operation .rot.

c = a .rot. b (where a, b, c are vec2d objects) :

• c%x = b%x * a%x + b%y * a%y
• c%y = b%x * a%y - b%y * a%x

The .invrot. operator points to the invrot_vec2d procedure. This operation unapplies the previously defined rotation operation.

Operation .invrot.

c = a .invrot. b (where a, b, c are vec2d objects) :

• c%x = b%x * a%x - b%y * a%y
• c%y = b%x * a%y + b%y * a%x

3D

The structure of a 3D vector contains three components x, y and z and several vector procedures. This structure overloads the assignment operator and other operators by creating a pointer component for each of them.

!===================================================================================================================!
!  3D Vector Structure
!===================================================================================================================!

TYPE vec3d

  real(rp) :: x = 0._rp
  real(rp) :: y = 0._rp
  real(rp) :: z = 0._rp

CONTAINS

  procedure, pass( self ) :: equal_vec3d
  procedure, pass( self ) :: equal_vec3d_scal
  procedure               :: add_vec3d
  procedure               :: sub_vec3d
  procedure, pass( self ) :: sub_vec3d_unary
  procedure, pass( self ) :: scalar_vec3d
  procedure, pass( self ) :: vec3d_scalar
  procedure               :: vec3d_div_scalar
  procedure               :: norm_vec3d
  procedure               :: square_vec3d
  procedure               :: dot_product_vec3d
  procedure               :: tens_prod3d
  procedure               :: rot_vec3d
  procedure               :: invrot_vec3d

  generic :: assignment(=) => equal_vec3d , equal_vec3d_scal

  generic :: operator(+) => add_vec3d
  generic :: operator(-) => sub_vec3d , sub_vec3d_unary
  generic :: operator(*) => scalar_vec3d , vec3d_scalar
  generic :: operator(/) => vec3d_div_scalar

  generic :: operator(.norm.    ) => norm_vec3d
  generic :: operator(.square.  ) => square_vec3d
  generic :: operator(.dotprod. ) => dot_product_vec3d
  generic :: operator(.tensprod.) => tens_prod3d
  generic :: operator(.rot.     ) => rot_vec3d
  generic :: operator(.invrot.  ) => invrot_vec3d

END TYPE vec3d

The assignement operator = points to the equal_vec3d and equal_vec3d_scal procedures.

The + operator points to the add_vec3d procedure.

The - operator points to the sub_vec3d and sub_vec3d_unary procedures.

The * operator points to the scalar_vec3d and vec3d_scalar procedures.

The / operator points to the vec3d_div_scalar procedure.

The .norm. operator points to the norm_vec3d procedure.

The .square. operator points to the square_vec3d procedure.

The .dotprod. operator points to the dot_product_vec3d procedure.

Assignment = and +, -, *, /, .norm., .square., .dotprod. operators

These operators are similar to the 2D operators, see Vectors - 2D to get the detailed operations and significations.

The .tensprod. operator points to the tens_prod3d procedure. This operation computes the tensor product of two vectors.

Operation .tensprod.

The .tensprod. operator is pretty similar to the 2D operator, still for c = a .tensprod. b (where a, b are vec3d objects and c a tens3d object) :

• c%xx = a%x * b%x
• c%xy = a%x * b%y
• c%xz = a%x * b%z
• c%yx = a%y * b%x
• c%yy = a%y * b%y
• c%yz = a%y * b%z
• c%zx = a%z * b%x
• c%zy = a%z * b%y
• c%zz = a%z * b%z

The .rot. operator points to the rot_vec3d procedure. This operation applies a rotation operation on a vector.

Operation .rot.

c = a .rot. b (where a, b, c are vec3d objects) :

d = sqrt( b%x**2 + b%y**2 )

if ( d < epsilon(1._rp) ) then

  c%x =   a%z
  c%y =   a%y
  c%z = - a%x

else

  dinv = - 1._rp / d

  c%x = b%x              * a%x + b%y              * a%y + b%z * a%z
  c%y = b%y       * dinv * a%x - b%x       * dinv * a%y
  c%z = b%x * b%z * dinv * a%x + b%y * b%z * dinv * a%y + d   * a%z

end if

The .invrot. operator points to the invrot_vec3d procedure. This operation unapplies the previously defined rotation operation.

Operation .invrot.

c = a .invrot. b (where a, b, c are vec3d objects) :

d = sqrt( b%x**2 + b%y**2 )

if ( d < epsilon(1._rp) ) then

  c%x = - a%z
  c%y =   a%y
  c%z =   a%x

else

  dinv = - 1._rp / d

  c%x = b%x * a%x + b%y * dinv * a%y + b%x * b%z * dinv * a%z
  c%y = b%y * a%x - b%x * dinv * a%y + b%y * b%z * dinv * a%z
  c%z = b%z * a%x                    + d                * a%z

end if

Tensor

Two structures for 2D and 3D tensors are implemented. The implementation of various operations is enhanced. One can easily add and substract tensors, multiply a tensor by a scalar or a vector, divide a tensor by a scalar, and transpose a tensor.

2D

The structure of a 2D tensor contains four components xx, xy, yx and yy and several tensor procedures. This structure overloads the assignment operator and other operators by creating a pointer component for each of them.

!===================================================================================================================!
!  2D Tensor Structure
!===================================================================================================================!

TYPE tens2d

  real(rp) :: xx = 0._rp
  real(rp) :: xy = 0._rp
  real(rp) :: yx = 0._rp
  real(rp) :: yy = 0._rp

CONTAINS

  procedure, pass( self ) :: equal_tens2d
  procedure, pass( self ) :: equal_tens2d_scal
  procedure               :: add_tens2d
  procedure               :: sub_tens2d
  procedure, pass( self ) :: scalar_tens2d
  procedure, pass( self ) :: tens2d_scalar
  procedure               :: tens2d_vec2d
  procedure               :: tens2d_div_scalar
  procedure               :: transpose_tens2d

  generic :: assignment(=) => equal_tens2d , equal_tens2d_scal

  generic :: operator(+) => add_tens2d
  generic :: operator(-) => sub_tens2d
  generic :: operator(*) => scalar_tens2d , tens2d_scalar , tens2d_vec2d
  generic :: operator(/) => tens2d_div_scalar

  generic :: operator(.t.) => transpose_tens2d

END TYPE tens2d

The assignement operator = points to the equal_tens2d and equal_tens2d_scal procedures. This operation assigns a value to each tensor component.

Assignement =

• a = b (where a and b are two tens2d objects) assigns all the b components to the corresponding a components.
• a = r (where a is a tens2d object and r is a scalar) assigns r to all the a components.

The + operator points to the add_tens2d procedure. This operation sums the components of each tensor.

Operation +

c = a + b (where a, b, c are tens2dobjects) :

• c%xx = a%xx + b%xx
• c%xy = a%xy + b%xy
• c%yx = a%yx + b%yx
• c%yy = a%yy + b%yy

The - operator points to the sub_tens2d procedure. This operation substract each component of one tensor to another.

Operation -

c = a - b (where a, b, c are tens2dobjects) :

• c%xx = a%xx - b%xx
• c%xy = a%xy - b%xy
• c%yx = a%yx - b%yx
• c%yy = a%yy - b%yy

The * operator points to the scalar_tens2d, tens2d_scalar, tens2d_vec2d procedures. This operation multiplies each tensor’s component by a scalar, or multiplies a tensor to a vector.

Operation *

• c = r * a or c = a * r (where a and c are tens2d objects, and r a scalar) :
  • c%xx = r * a%xx
  • c%xy = r * a%xy
  • c%yx = r * a%yx
  • c%yy = r * a%yy
• c = A * b (where A is a tens2d object and c and b vec2dobjects) :
  • c%x = A%xx * b%x + A%xy * b%y
  • c%y = A%yx * b%x + A%yy * b%y

The / operator points to the tens2d_div_scalar procedure. This operation divides each tensor’s component by a scalar.

Operation /

c = a / r (where a and c are tens2d objects and r a scalar) :

• c%xx = a%xx / r
• c%xy = a%xy / r
• c%yx = a%yx / r
• c%yy = a%yy / r

The .t. operator points to the transpose_tens2d procedure. This operation computes the tensor’s transposition.

Operation .t.

c = .t. a (where a and c are tens2d objects) :

• c%xx = a%xx
• c%xy = a%yx
• c%yx = a%xy
• c%yy = a%yy

3D

The structure of a 3D tensor contains nine components xx, xy, xz, yx, yy, yz, zx, zy and zz and several tensor procedures. This structure overloads the assignment operator and other operators by creating a pointer component for each of them.

!===================================================================================================================!
!  3D Tensor Structure
!===================================================================================================================!

TYPE tens3d

  real(rp) :: xx = 0._rp
  real(rp) :: xy = 0._rp
  real(rp) :: xz = 0._rp
  real(rp) :: yx = 0._rp
  real(rp) :: yy = 0._rp
  real(rp) :: yz = 0._rp
  real(rp) :: zx = 0._rp
  real(rp) :: zy = 0._rp
  real(rp) :: zz = 0._rp

CONTAINS

  procedure, pass( self ) :: equal_tens3d
  procedure, pass( self ) :: equal_tens3d_scal
  procedure               :: add_tens3d
  procedure               :: sub_tens3d
  procedure, pass( self ) :: scalar_tens3d
  procedure, pass( self ) :: tens3d_scalar
  procedure               :: tens3d_vec3d
  procedure               :: tens3d_div_scalar
  procedure               :: transpose_tens3d

  generic :: assignment(=) => equal_tens3d , equal_tens3d_scal

  generic :: operator(+) => add_tens3d
  generic :: operator(-) => sub_tens3d
  generic :: operator(*) => scalar_tens3d , tens3d_scalar , tens3d_vec3d
  generic :: operator(/) => tens3d_div_scalar

  generic :: operator(.t.) => transpose_tens3d

END TYPE tens3d

The assignement operator = points to the equal_tens3d and equal_tens3d_scal procedures.

The + operator points to the add_tens3d procedure.

The - operator points to the sub_tens3d procedure.

Assignment = and +, - operators

These operators are similar to the 2D operators, see Tensors - 2D to get the detailed operations and significations.

The * operator points to the scalar_tens3d, tens3d_scalar, tens3d_vec3d procedures. This operation multiplies each tensor’s component by a scalar, or multiplies a tensor to a vector.

Operation *

• c = r * a or c = a * r : If a scalar is involved in the operation, the multiplication is similar to the 2D multiplication (see Tensors - 2D).

• c = A * b (where A is a tens3d object and c and b vec3dobjects) :

  • c%x = A%xx * b%x + A%xy * b%y + A%xz * b%z
  • c%y = A%yx * b%x + A%yy * b%y + A%yz * b%z
  • c%z = A%zx * b%x + A%zy * b%y + A%zz * b%z

The / operator points to the tens3d_div_scalar procedure.

Operation /

This operator is similar to the 2D operator, see Tensors - 2D to get the detailed operation and signification.

The .t. operator points to the transpose_tens3d procedure. This operation computes the tensor’s transposition.

Operation .t.

c = .t. a (where a and c are tens3d objects) :

• c%xx = a%xx
• c%xy = a%yx
• c%xz = a%zx
• c%yx = a%xy
• c%yy = a%yy
• c%yz = a%zy
• c%zx = a%xz
• c%zy = a%yz
• c%zz = a%zz