mixedpoisson.h

Go to the documentation of this file.

//
//  mixedpoisson.h
//  PZ
//
//  Created by Agnaldo Farias on 5/28/12.
//  Copyright (c) 2012 LabMec-Unicamp. All rights reserved.
//

#ifndef PZ_mixedpoisson_h
#define PZ_mixedpoisson_h

#include "TPZMaterial.h"
#include "pzdiscgal.h"
#include "pzpoisson3d.h"
#include "TPZMaterial.h"
#include "pzfunction.h"

class TPZMixedPoisson : public TPZMatPoisson3d {
    
protected:
                REAL ff;
                
                TPZFNMatrix<9,REAL> fTensorK;
    
                TPZFNMatrix<9,REAL> fInvK;
    
                REAL fvisc;
    
    bool fIsStabilized;
    
    REAL fdelta1;
    REAL fdelta2;
    
    REAL fh2;
    bool fUseHdois;
    TPZAutoPointer<TPZFunction<STATE> > fPermeabilityFunction;
    
public:
    TPZMixedPoisson();
    
    TPZMixedPoisson(int matid, int dim);
    
    virtual ~TPZMixedPoisson();
    
    TPZMixedPoisson(const TPZMixedPoisson &cp);
    
    TPZMixedPoisson &operator=(const TPZMixedPoisson &copy);
    
    virtual TPZMaterial * NewMaterial() override{
        return new TPZMixedPoisson(*this);
    }
    

    virtual void Print(std::ostream & out) override;
                
                virtual std::string Name() override{ return "TPZMixedPoisson"; }
    
    virtual int NStateVariables() const override;
                
                void SetPermeability(REAL perm) {
                                fK = perm;
        fTensorK.Zero();
        fInvK.Zero();
        for (int i=0; i<3; i++) {
            fTensorK(i,i) = perm;
            fInvK(i,i) = 1./perm;
        }
                }
    
    //Set the permeability tensor and inverser tensor
   void SetPermeabilityTensor(const TPZFMatrix<REAL> &K, const TPZFMatrix<REAL> &invK){
   
//       if(K.Rows() != fDim || K.Cols() != fDim) DebugStop();
//       if(K.Rows()!=invK.Rows() || K.Cols()!=invK.Cols()) DebugStop();
       if(K.Rows() != 3 || invK.Rows() != 3)
       {
           DebugStop();
       }
       fTensorK = K;
       fInvK = invK;
       fK = 0.;
   }
                
    void GetPermeability(TPZVec<REAL> &x, TPZFMatrix<REAL> &K, TPZFMatrix<REAL> &invK);

    void SetViscosity(REAL visc) {
                                fvisc = visc;
                }
    
    void GetMaxPermeability(REAL &perm)
    {
        perm = 0;
        for(int i=0; i<3; i++) perm = perm < fTensorK(i,i) ? fTensorK(i,i) : perm;
    }
//              void GetPermeability(REAL &perm) {
//                              perm = fK;
//              }
                
                void SetInternalFlux(REAL flux) {
                                ff = flux;
                }
    
    void SetStabilizedMethod(){
        fIsStabilized = true;
    }
    
    void SetHdois(){
        fUseHdois = true;
    }
    
    void SetStabilizationCoeficients(REAL delta1, REAL delta2){
        fdelta1 = delta1;
        fdelta2 = delta2;
    }
    
    void SetPermeabilityFunction(TPZAutoPointer<TPZFunction<STATE> > fp)
    {
        fPermeabilityFunction = fp;
    }

    TPZAutoPointer<TPZFunction<STATE> > PermeabilityFunction()
    {
        return fPermeabilityFunction;
    }
    
    virtual void Contribute(TPZVec<TPZMaterialData> &datavec, REAL weight, TPZFMatrix<STATE> &ek, TPZFMatrix<STATE> &ef) override;
                
    virtual void ContributeBC(TPZVec<TPZMaterialData> &datavec, REAL weight, TPZFMatrix<STATE> &ek,TPZFMatrix<STATE> &ef,TPZBndCond &bc) override;

                
                virtual int VariableIndex(const std::string &name) override;
                
                virtual int NSolutionVariables(int var) override;
                
     virtual void Solution(TPZVec<TPZMaterialData> &datavec, int var, TPZVec<STATE> &Solout) override;
    
    virtual void FillBoundaryConditionDataRequirement(int type,TPZVec<TPZMaterialData > &datavec) override
    {
        // default is no specific data requirements
        int nref = datavec.size();
        for (int iref = 0; iref <nref; iref++) {
            datavec[iref].SetAllRequirements(false);
            datavec[iref].fNeedsSol = true;
        }
        datavec[0].fNeedsNormal = true;
        if(type == 50)
        {
            for(int iref = 0; iref<nref; iref++){
                datavec[iref].fNeedsSol = true;
            }
        }
    }

    virtual void FillDataRequirements(TPZVec<TPZMaterialData > &datavec) override;
                
    
    virtual int NEvalErrors() override {return 3;}

    virtual void Errors(TPZVec<TPZMaterialData> &data, TPZVec<STATE> &u_exact, TPZFMatrix<STATE> &du_exact, TPZVec<REAL> &errors) override;
    
    public:
virtual int ClassId() const  override;

};

#endif