huertaerrorestimator.h

Go to the documentation of this file.

/*
 *
 *                 #####    #####   ######  ######  ###   ###
 *               ##   ##  ##   ##  ##      ##      ## ### ##
 *              ##   ##  ##   ##  ####    ####    ##  #  ##
 *             ##   ##  ##   ##  ##      ##      ##     ##
 *            ##   ##  ##   ##  ##      ##      ##     ##
 *            #####    #####   ##      ######  ##     ##
 *
 *
 *             OOFEM : Object Oriented Finite Element Code
 *
 *               Copyright (C) 1993 - 2013   Borek Patzak
 *
 *
 *
 *       Czech Technical University, Faculty of Civil Engineering,
 *   Department of Structural Mechanics, 166 29 Prague, Czech Republic
 *
 *  This library is free software; you can redistribute it and/or
 *  modify it under the terms of the GNU Lesser General Public
 *  License as published by the Free Software Foundation; either
 *  version 2.1 of the License, or (at your option) any later version.
 *
 *  This program 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
 *  Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public
 *  License along with this library; if not, write to the Free Software
 *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
 */

#ifndef huertaerrorestimator_h
#define huertaerrorestimator_h

#include "../sm/refinedelement.h"
#include "../sm/refinedmesh.h"
#include "errorestimator.h"
#include "interface.h"
#include "floatarray.h"
#include "statecountertype.h"
#include "remeshingcrit.h"

#include <vector>


#define _IFT_HuertaErrorEstimator_Name "huerta"
#define _IFT_HuertaErrorEstimator_normtype "normtype"
#define _IFT_HuertaErrorEstimator_refinelevel "refinelevel"
#define _IFT_HuertaErrorEstimator_requirederror "requirederror"
#define _IFT_HuertaErrorEstimator_skipsteps "skipsteps"
#define _IFT_HuertaErrorEstimator_initialskipsteps "initialskipsteps"
#define _IFT_HuertaErrorEstimator_werror "werror"
#define _IFT_HuertaErrorEstimator_perfectCSect "perCSect"
#define _IFT_HuertaErrorEstimator_impCSect "impCSect"
#define _IFT_HuertaErrorEstimator_impPos "imppos"
#define _IFT_HuertaErrorEstimator_exact "exact"



#define _IFT_HuertaRemeshingCriteria_requirederror "requirederror"
#define _IFT_HuertaRemeshingCriteria_minelemsize "minelemsize"
#define _IFT_HuertaRemeshingCriteria_noremesh "noremesh"
#define _IFT_HuertaRemeshingCriteria_werror "werror"
#define _IFT_HuertaRemeshingCriteria_refinecoeff "refinecoeff"



#define _IFT_HuertaErrorEstimatorInterface_coords "coords"


namespace oofem {
class Element;
class GaussPoint;

class HuertaErrorEstimator : public ErrorEstimator
{
public:
    enum NormType { L2Norm, EnergyNorm };
    enum AnalysisMode { HEE_linear, HEE_nlinear };

protected:
    double globalENorm;
    double globalWENorm;
    double globalUNorm;
    double globalErrorEstimate;
    FloatArray eNorms;
    NormType normType;
    StateCounterType stateCounter;
    FloatArray primaryUnknownError;
    int refineLevel;
    std :: vector< RefinedElement >refinedElementList;
    RefinedMesh refinedMesh;
    AnalysisMode mode;
    double requiredError;
    bool wError;

    double lastError;
    int stepsToSkip, skippedSteps, maxSkipSteps, initialSkipSteps;

public:
    HuertaErrorEstimator(int n, Domain * d) : ErrorEstimator(n, d), eNorms(0), primaryUnknownError(0),
        refinedElementList(), refinedMesh()
    {
        eeType = EET_HEE;
        stateCounter = 0;
        normType = EnergyNorm;
        refineLevel = 1;
        mode = HEE_linear;
        wError = false;
        lastError = -1.0;
        stepsToSkip = skippedSteps = initialSkipSteps = 0;
    }

    virtual ~HuertaErrorEstimator() { }

    int giveRefinementLevel() { return this->refineLevel; }

    virtual double giveElementError(EE_ErrorType type, Element *elem, TimeStep *tStep);

    virtual double giveValue(EE_ValueType type, TimeStep *tStep);

    virtual int estimateError(EE_ErrorMode err_mode, TimeStep *tStep);

    virtual RemeshingCriteria *giveRemeshingCrit();

    virtual IRResultType initializeFrom(InputRecord *ir);

    virtual const char *giveInputRecordName() const { return NULL; }
    virtual const char *giveClassName() const { return "HuertaErrorEstimator"; }

    AnalysisMode giveAnalysisMode() { return mode; }

    virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj = NULL);
    virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj = NULL);

private:
    void buildRefinedMesh();

    void solveRefinedElementProblem(int elemId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
                                    TimeStep *tStep);
    void solveRefinedPatchProblem(int nodeId, IntArray &localNodeIdArray,
                                  IntArray &globalNodeIdArray, TimeStep *tStep);
    void solveRefinedWholeProblem(IntArray &localNodeIdArray, IntArray &globalNodeIdArray, TimeStep *tStep);
    void extractVectorFrom(Element *element, FloatArray &vector, FloatArray &answer, int dofs, TimeStep *tStep);

    void setupRefinedProblemProlog(const char *problemName, int problemId, IntArray &localNodeIdArray,
                                   int nodes, int elems, int csects, int mats, int loads, int funcs,
                                   IntArray &controlNode, IntArray &controlDof, TimeStep *tStep);
    void setupRefinedProblemEpilog1(int csects, int mats, int loads, int nlbarriers);
    void setupRefinedProblemEpilog2(int tfuncs);
};


class HuertaErrorEstimatorInterface : public Interface
{
public:
    enum SetupMode { CountMode = 0, NodeMode = 1, ElemMode = 2, BCMode = 3 };

public:
    HuertaErrorEstimatorInterface() { }

    virtual void HuertaErrorEstimatorI_setupRefinedElementProblem(RefinedElement *refinedElement, int level, int nodeId,
                                                                  IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
                                                                  HuertaErrorEstimatorInterface :: SetupMode mode, TimeStep *tStep,
                                                                  int &localNodeId, int &localElemId, int &localBcId,
                                                                  IntArray &controlNode, IntArray &controlDof,
                                                                  HuertaErrorEstimator :: AnalysisMode aMode) = 0;

    virtual void HuertaErrorEstimatorI_computeNmatrixAt(GaussPoint *gp, FloatMatrix &answer) = 0;

protected:
    void setupRefinedElementProblem1D(Element *element, RefinedElement *refinedElement,
                                      int level, int nodeId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
                                      HuertaErrorEstimatorInterface :: SetupMode mode, TimeStep *tStep, int nodes,
                                      FloatArray **corner, FloatArray &midNode,
                                      int &localNodeId, int &localElemId, int &localBcId,
                                      IntArray &controlNode, IntArray &controlDof,
                                      HuertaErrorEstimator :: AnalysisMode aMode, const char *edgetype);

    void setupRefinedElementProblem2D(Element *element, RefinedElement *refinedElement,
                                      int level, int nodeId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
                                      HuertaErrorEstimatorInterface :: SetupMode mode, TimeStep *tStep, int nodes,
                                      FloatArray **corner, FloatArray *midSide, FloatArray &midNode,
                                      int &localNodeId, int &localElemId, int &localBcId,
                                      IntArray &controlNode, IntArray &controlDof,
                                      HuertaErrorEstimator :: AnalysisMode aMode, const char *quadtype);

    void setupRefinedElementProblem3D(Element *element, RefinedElement *refinedElement,
                                      int level, int nodeId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
                                      HuertaErrorEstimatorInterface :: SetupMode mode, TimeStep *tStep, int nodes,
                                      FloatArray **corner, FloatArray *midSide, FloatArray *midFace, FloatArray &midNode,
                                      int &localNodeId, int &localElemId, int &localBcId,
                                      int hexaSideNode [ 1 ] [ 3 ], int hexaFaceNode [ 1 ] [ 3 ],
                                      IntArray &controlNode, IntArray &controlDof,
                                      HuertaErrorEstimator :: AnalysisMode aMode, const char *hexatype);
};


class HuertaRemeshingCriteria : public RemeshingCriteria
{
public:
    enum HuertaRemeshingCriteriaModeType { primaryUnknownBased };

protected:
    FloatArray nodalDensities;
    RemeshingStrategy remeshingStrategy;
    StateCounterType stateCounter;
    HuertaRemeshingCriteriaModeType mode;
    double requiredError;
    double minElemSize;
    double refineCoeff;
    bool noRemesh;
    bool wError;

public:
    HuertaRemeshingCriteria(int n, ErrorEstimator * e);
    virtual ~HuertaRemeshingCriteria() { }

    virtual double giveRequiredDofManDensity(int num, TimeStep *tStep, int relative = 0);
    virtual double giveDofManDensity(int num);
    virtual RemeshingStrategy giveRemeshingStrategy(TimeStep *tStep);
    virtual int estimateMeshDensities(TimeStep *tStep);
    virtual IRResultType initializeFrom(InputRecord *ir);

    virtual const char *giveInputRecordName() const { return NULL; }
    virtual const char *giveClassName() const { return "HuertaErrorEstimator"; }
};

} // end namespace oofem
#endif // huertaerrorestimator_h