geometry.h

Go to the documentation of this file.

#ifndef GEOMETRY_H
#define GEOMETRY_H
 
#include "apps/wmapplication.h" // for solverVariables_t
 
inline void rotate_vector(const real v[3], const solverVariables_t* pSV, real Rv[3])
{
 
    Rv[0] = pSV->R[0][0] * v[0] + pSV->R[0][1] * v[1] + pSV->R[0][2] * v[2];
    Rv[1] = pSV->R[1][0] * v[0] + pSV->R[1][1] * v[1] + pSV->R[1][2] * v[2];
    Rv[2] = pSV->R[2][0] * v[0] + pSV->R[2][1] * v[1] + pSV->R[2][2] * v[2];
}
 
inline void antirotate_vector(const real Rv[3], const solverVariables_t* pSV, real v[3])
{
    v[0] = pSV->R[0][0] * Rv[0] + pSV->R[1][0] * Rv[1] + pSV->R[2][0] * Rv[2];
    v[1] = pSV->R[0][1] * Rv[0] + pSV->R[1][1] * Rv[1] + pSV->R[2][1] * Rv[2];
    v[2] = pSV->R[0][2] * Rv[0] + pSV->R[1][2] * Rv[1] + pSV->R[2][2] * Rv[2];
}
 
inline void rotate_matrix(const real* m[3], const solverVariables_t* pSV, real* RmRT[3])
{
    // RmRT = R \cdot m \cdot R^T
    real val;
    for (int i = 0; i < 3; i++)
    {
        for (int j = 0; j < 3; j++)
        {
            val = 0.;
            for (int k = 0; k < 3; k++)
            {
                for (int l = 0; l < 3; l++)
                {
                    val += pSV->R[i][k] * m[k][l] * pSV->R[j][l];
                }
            }
            RmRT[i][j] = val;
        }
    }
}
 
inline void
antirotate_matrix(const real* RmRT[3], const solverVariables_t* pSV, real* m[3])
{
    // m = R^T \cdot RmRT \cdot R
    real val;
    for (int i = 0; i < 3; i++)
    {
        for (int j = 0; j < 3; j++)
        {
            val = 0.;
            for (int k = 0; k < 3; k++)
            {
                for (int l = 0; l < 3; l++)
                {
                    val += pSV->R[k][i] * RmRT[k][l] * pSV->R[l][j];
                }
            }
            m[i][j] = val;
        }
    }
}
 
inline real dot_product(const real* v)
{
    return (v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
}
 
inline real dot_product(const real* u, const real* v)
{
    return (u[0] * v[0] + u[1] * v[1] + u[2] * v[2]);
}
 
inline real dot_product(const real x, const real y, const real z)
{
    return (x * x + y * y + z * z);
}
 
inline real dot_product(const real x1,
                        const real y1,
                        const real z1,
                        const real x2,
                        const real y2,
                        const real z2)
{
    return (x1 * x2 + y1 * y2 + z1 * z2);
}
 
inline real levi_civita(size_t a, size_t b, size_t c)
{
    if (a == b || b == c || c == a)
    {
        return 0;
    }
    else if (b == (a + 1) || c == (b + 1) || a == (c + 1))
    {
        return 1;
    }
    return -1;
}
 
#endif // GEOMETRY_H