/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | www.openfoam.com
     \\/     M anipulation  |
-------------------------------------------------------------------------------
    Copyright (C) 2016 OpenFOAM Foundation
    Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Class
    Foam::MatrixSpace

Description
    Templated matrix space.

    Template arguments are the Form the matrix space will be used to create,
    the type of the elements and the number of rows and columns of the matrix.

SourceFiles
    MatrixSpaceI.H

See also
    Foam::VectorSpace

\*---------------------------------------------------------------------------*/

#ifndef MatrixSpace_H
#define MatrixSpace_H

#include "VectorSpace.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

/*---------------------------------------------------------------------------*\
                         Class MatrixSpace Declaration
\*---------------------------------------------------------------------------*/

template<class Form, class Cmpt, direction Mrows, direction Ncols>
class MatrixSpace
:
    public VectorSpace<Form, Cmpt, Mrows*Ncols>
{
public:

    // Typedefs

        //- MatrixSpace type
        typedef MatrixSpace<Form, Cmpt, Mrows, Ncols> msType;


    // Member Constants

        static constexpr direction mRows = Mrows;
        static constexpr direction nCols = Ncols;


    // Static Member Functions

        //- The number of rows
        __host__ __device__
        static direction m() noexcept
        {
            return Mrows;
        }

        //- The number of columns
        __host__ __device__
        static direction n() noexcept
        {
            return Ncols;
        }

        //- An identity matrix for square matrix-spaces
        __host__ __device__
        inline static msType identity();


    // Sub-Block Classes

        //- Const sub-block type
        template<class SubTensor, direction BRowStart, direction BColStart>
        class ConstBlock
        {
            //- Reference to parent matrix
            const msType& matrix_;

        public:

            static const direction mRows = SubTensor::mRows;
            static const direction nCols = SubTensor::nCols;

            //- Return the number of rows in the block
            __host__ __device__
            static direction m()
            {
                return mRows;
            }

            //- Return the number of columns in the block
            __host__ __device__
            static direction n()
            {
                return nCols;
            }

            //- Construct for the given matrix
            __host__ __device__
            inline ConstBlock(const msType& matrix);

            //- Construct and return the sub-tensor corresponding to this block
            __host__ __device__
            inline SubTensor operator()() const;

            //- (i, j) const element access operator
            __host__ __device__
            inline const Cmpt& operator()
            (
                const direction i,
                const direction j
            ) const;
        };


        //- Sub-block type
        template
        <
            class SubTensor,
            direction BRowStart,
            direction BColStart
        >
        class Block
        {
            //- Reference to parent matrix
            msType& matrix_;

        public:

            static const direction mRows = SubTensor::mRows;
            static const direction nCols = SubTensor::nCols;

            //- The number of rows in the block
            __host__ __device__
            static direction m()
            {
                return mRows;
            }

            //- The number of columns in the block
            __host__ __device__
            static direction n()
            {
                return nCols;
            }

            //- Construct for the given matrix
            __host__ __device__
            inline Block(msType& matrix);

            //- Assignment to a matrix
            template<class Form2>
            __host__ __device__
            inline void operator=
            (
                const MatrixSpace
                <
                    Form2,
                    Cmpt,
                    SubTensor::mRows,
                    SubTensor::nCols
                >& matrix
            );

            //- Assignment to a column vector
            template<class VSForm>
            __host__ __device__
            inline void operator=
            (
                const VectorSpace<VSForm, Cmpt, SubTensor::mRows>& v
            );

            //- Construct and return the sub-tensor corresponding to this block
            __host__ __device__
            inline SubTensor operator()() const;

            //- (i, j) const element access operator
            __host__ __device__
            inline const Cmpt& operator()
            (
                const direction i,
                const direction j
            ) const;

            //- (i, j) element access operator
            __host__ __device__
            inline Cmpt& operator()(const direction i, const direction j);
        };


    // Generated Methods

        //- Default construct
        __host__ __device__
        MatrixSpace() = default;


    // Constructors

        //- Construct initialized to zero
        __host__ __device__
        inline MatrixSpace(const Foam::zero);

        //- Construct as copy of a VectorSpace with the same size
        template<class Form2, class Cmpt2>
        __host__ __device__
        inline explicit MatrixSpace
        (
            const VectorSpace<Form2, Cmpt2, Mrows*Ncols>&
        );

        //- Construct from a block of another matrix space
        template
        <
            template<class, direction, direction> class Block2,
            direction BRowStart,
            direction BColStart
        >
        __host__ __device__
        inline MatrixSpace
        (
            const Block2<Form, BRowStart, BColStart>& block
        );

        //- Construct from Istream
        __host__ __device__
        explicit MatrixSpace(Istream& is);


    // Member Functions

        //- Fast const element access using compile-time addressing
        template<direction Row, direction Col>
        __host__ __device__
        inline const Cmpt& elmt() const;

        //- Fast element access using compile-time addressing
        template<direction Row, direction Col>
        __host__ __device__ 
        inline Cmpt& elmt();

        // Const element access functions for a 3x3
        // Compile-time errors are generated for inappropriate use

            __host__ __device__
            inline const Cmpt& xx() const;
            __host__ __device__
            inline const Cmpt& xy() const;
            __host__ __device__
            inline const Cmpt& xz() const;
            __host__ __device__
            inline const Cmpt& yx() const;
            __host__ __device__
            inline const Cmpt& yy() const;
            __host__ __device__
            inline const Cmpt& yz() const;
            __host__ __device__
            inline const Cmpt& zx() const;
            __host__ __device__
            inline const Cmpt& zy() const;
            __host__ __device__
            inline const Cmpt& zz() const;

        // Element access functions for a 3x3
        // Compile-time errors are generated for inappropriate use

            __host__ __device__
            inline Cmpt& xx();
            __host__ __device__
            inline Cmpt& xy();
            __host__ __device__
            inline Cmpt& xz();
            __host__ __device__
            inline Cmpt& yx();
            __host__ __device__
            inline Cmpt& yy();
            __host__ __device__
            inline Cmpt& yz();
            __host__ __device__
            inline Cmpt& zx();
            __host__ __device__
            inline Cmpt& zy();
            __host__ __device__
            inline Cmpt& zz();

        //- Return the transpose of the matrix
        __host__ __device__
        inline typename typeOfTranspose<Cmpt, Form>::type T() const;

        //- Return a const sub-block corresponding to the specified type
        //  starting at the specified row and column
        template<class SubTensor, direction BRowStart, direction BColStart>
        __host__ __device__
        inline ConstBlock<SubTensor, BRowStart, BColStart> block() const;

        //- Return a sub-block corresponding to the specified type
        //  starting at the specified row and column
        template<class SubTensor, direction BRowStart, direction BColStart>
        __host__ __device__
        inline Block<SubTensor, BRowStart, BColStart> block();

        //- (i, j) const element access operator
        __host__ __device__
        inline const Cmpt& operator()
        (
            const direction& i,
            const direction& j
        ) const;

        //- (i, j) element access operator
        __host__ __device__
        inline Cmpt& operator()(const direction& i, const direction& j);


    // Member Operators

        //- Assignment to zero
        __host__ __device__
        inline void operator=(const Foam::zero);

        //- Assignment to a block of another matrix space
        template
        <
            template<class, direction, direction> class Block2,
            direction BRowStart,
            direction BColStart
        >
        __host__ __device__
        inline void operator=
        (
            const Block2<Form, BRowStart, BColStart>& block
        );

        //- Inner product with a compatible square matrix
        template<class Form2>
        __host__ __device__
        inline void operator&=
        (
            const MatrixSpace<Form, Cmpt, Ncols, Ncols>& matrix
        );
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#include "MatrixSpaceI.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //
