MultisegmentWellContribution.hpp

/*
  Copyright 2020 Equinor ASA

  This file is part of the Open Porous Media project (OPM).

  OPM 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.

  OPM 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 OPM.  If not, see <http://www.gnu.org/licenses/>.
*/

#ifndef MULTISEGMENTWELLCONTRIBUTION_HEADER_INCLUDED
#define MULTISEGMENTWELLCONTRIBUTION_HEADER_INCLUDED

#include <vector>

#if HAVE_CUDA
#include <cuda_runtime.h>
#endif

#include <umfpack.h>
#include <dune/common/version.hh>

namespace Opm
{


template<class Scalar>
class MultisegmentWellContribution
{

private:
    unsigned int dim;                        // size of blockvectors in vectors x and y, equal to MultisegmentWell::numEq
    unsigned int dim_wells;                  // size of blocks in C, B and D, equal to MultisegmentWell::numWellEq
    unsigned int M;                          // number of rows, M == dim_wells*Mb
    unsigned int Mb;                         // number of blockrows in C, D and B

#if HAVE_CUDA
    cudaStream_t stream; // not actually used yet, will be when MultisegmentWellContribution are applied on GPU
#endif

    // C and B are stored in BCRS format, D is stored in CSC format (Dune::UMFPack)
    // Sparsity pattern for C is not stored, since it is the same as B
    unsigned int DnumBlocks;             // number of blocks in D
    std::vector<Scalar> Cvals;
    std::vector<Scalar> Dvals;
    std::vector<Scalar> Bvals;
    std::vector<int> Dcols;              // Columnpointers, contains M+1 entries
    std::vector<unsigned int> Bcols;
    std::vector<int> Drows;              // Rowindicies, contains DnumBlocks*dim*dim_wells entries
    std::vector<unsigned int> Brows;
    std::vector<Scalar> z1;          // z1 = B * x
    std::vector<Scalar> z2;          // z2 = D^-1 * B * x
    void *UMFPACK_Symbolic, *UMFPACK_Numeric;

    unsigned int getColIdx(unsigned int idx);

public:
    using UMFPackIndex = SuiteSparse_long;

#if HAVE_CUDA
    void setCudaStream(cudaStream_t stream);
#endif

    MultisegmentWellContribution(unsigned int dim, unsigned int dim_wells,
                                 unsigned int Mb,
                                 std::vector<Scalar>& Bvalues,
                                 std::vector<unsigned int>& BcolIndices,
                                 std::vector<unsigned int>& BrowPointers,
                                 unsigned int DnumBlocks,
                                 Scalar* Dvalues,
                                 UMFPackIndex* DcolPointers,
                                 UMFPackIndex* DrowIndices,
                                 std::vector<Scalar>& Cvalues);

    ~MultisegmentWellContribution();

    void apply(Scalar* h_x, Scalar* h_y);
};

} //namespace Opm

#endif