source: issm/trunk/src/c/objects/FemModel.cpp@ 4573

/*!\file FemModel.c
 * \brief: implementation of the FemModel object
 */


#ifdef HAVE_CONFIG_H
        #include "config.h"
#else
#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
#endif

#include "stdio.h"
#include "../Container/Container.h"
#include "../modules/ModelProcessorx/ModelProcessorx.h"
#include "./objects.h"
#include "../include/include.h"
#include "../EnumDefinitions/EnumDefinitions.h"
#include "../modules/modules.h"

/*Object constructors and destructor*/
/*FUNCTION FemModel::constructor {{{1*/
FemModel::FemModel(ConstDataHandle IOMODEL,const int in_solution_type,const int* analyses,const int nummodels){

        /*intermediary*/
        int i;
        IoModel* iomodel=NULL;
        int analysis_type;

        /*Initialize internal data: */
        this->nummodels=nummodels;
        this->solution_type=in_solution_type;
        this->analysis_counter=nummodels-1; //point to last analysis_type carried out.
        this->results=new DataSet(); //not initialized by CreateDataSets
        this->partition=NULL;
        this->tpartition=NULL;
        
        /*Dynamically allocate whatever is a list of length nummodels: */
        analysis_type_list=(int*)xmalloc(nummodels*sizeof(int));
        m_Rmg=(Mat*)xmalloc(nummodels*sizeof(Mat));
        m_Gmn=(Mat*)xmalloc(nummodels*sizeof(Mat));
        m_nodesets=(NodeSets**)xmalloc(nummodels*sizeof(NodeSets*));
        m_yg=(Vec*)xmalloc(nummodels*sizeof(Vec));
        m_ys=(Vec*)xmalloc(nummodels*sizeof(Vec));

        /*Initialize: */
        for(i=0;i<nummodels;i++)analysis_type_list[i]=analyses[i];
        for(i=0;i<nummodels;i++)m_Rmg[i]=NULL;
        for(i=0;i<nummodels;i++)m_Gmn[i]=NULL;
        for(i=0;i<nummodels;i++)m_nodesets[i]=NULL;
        for(i=0;i<nummodels;i++)m_yg[i]=NULL;
        for(i=0;i<nummodels;i++)m_ys[i]=NULL;

        /*create datasets for all analyses: */
        ModelProcessorx(&this->elements,&this->nodes,&this->vertices,&this->materials,&this->constraints,&this->loads,&this->parameters,IOMODEL,this->solution_type,nummodels,analyses);

        /*do the post-processing of the datasets to get an FemModel that can actually run analyses: */
        for(i=0;i<nummodels;i++){

                _printf_("   processing finite element model of analysis %s:\n",EnumAsString(analysis_type_list[i]));
                analysis_type=analysis_type_list[i];
                this->SetCurrentConfiguration(analysis_type);
        
                _printf_("      create degrees of freedom\n");
                VerticesDofx(&partition,&tpartition,vertices,parameters);
                NodesDofx(nodes,parameters,analysis_type);

                _printf_("      create single point constraints\n");
                SpcNodesx( &m_yg[i], nodes,constraints,analysis_type); 

                _printf_("      create rigid body constraints\n");
                MpcNodesx( &m_Rmg[i], nodes,constraints,analysis_type); 

                _printf_("      create node sets\n");
                BuildNodeSetsx(&m_nodesets[i], nodes,analysis_type);

                _printf_("      reducing single point constraints vector\n");
                Reducevectorgtosx(&m_ys[i], m_yg[i],m_nodesets[i]);

                _printf_("      normalizing rigid body constraints matrix\n");
                NormalizeConstraintsx(&m_Gmn[i], m_Rmg[i],m_nodesets[i]);

                _printf_("      configuring element and loads\n");
                ConfigureObjectsx(elements, loads, nodes, vertices, materials,parameters);
        }
}
/*}}}1*/
/*FUNCTION FemModel::destructor {{{1*/
FemModel::~FemModel(){

        /*Intermediary*/
        int i;

        /*Delete all the datasets: */
        xfree((void**)&analysis_type_list);
        delete elements;
        delete nodes;
        delete vertices;
        delete constraints;
        delete loads;
        delete materials;
        delete parameters;
        delete results;
        VecFree(&partition);
        VecFree(&tpartition);

        for(i=0;i<nummodels;i++){
                Mat temp_Rmg=m_Rmg[i];
                MatFree(&temp_Rmg);
                Mat temp_Gmn=m_Gmn[i];
                MatFree(&temp_Gmn);
                NodeSets* temp_nodesets=m_nodesets[i];
                delete temp_nodesets;
                Vec temp_yg=m_yg[i];
                VecFree(&temp_yg);
                Vec temp_ys=m_ys[i];
                VecFree(&temp_ys);
        }

        /*Delete dynamically allocated arrays: */
        xfree((void**)&m_Rmg);
        xfree((void**)&m_Gmn);
        xfree((void**)&m_nodesets);
        xfree((void**)&m_yg);
        xfree((void**)&m_ys);

}
/*}}}1*/

/*Object management*/
/*FUNCTION FemModel::Echo {{{1*/
void FemModel::Echo(void){

        printf("FemModel echo: \n");
        printf("   number of fem models: %i\n",nummodels);
        printf("   analysis_type_list: \n");
        for(int i=0;i<nummodels;i++)printf("     %i: %s\n",i,EnumAsString(analysis_type_list[i]));
        printf("   current analysis_type: \n");
        printf("     %i: %s\n",analysis_counter,EnumAsString(analysis_type_list[analysis_counter]));


}
/*}}}*/

/*Numerics: */
/*FUNCTION FemModel::SetCurrentConfiguration(int configuration_type,int analysis_type){{{1*/
void FemModel::SetCurrentConfiguration(int configuration_type,int analysis_type){

        /*Use configuration_type to setup the analysis counter, the configurations of objects etc ... but use 
         * analysis_type to drive the element numerics. This allows for use of 1 configuration_type for several 
         * analyses. For example: do a SurfaceSlopeX, SurfaceSlopeY, BedSlopeX and BedSlopeY analysis using the 
         * Slope configuration.*/

        int found=-1;
        for(int i=0;i<nummodels;i++){
                if (analysis_type_list[i]==configuration_type){
                        found=i;
                        break;
                }
        }
        if(found!=-1) analysis_counter=found;
        else ISSMERROR("Could not find alias for analysis_type %s in list of FemModel analyses",EnumAsString(configuration_type));

        /*activate matrices/vectors: */
        Rmg=m_Rmg[analysis_counter];
        Gmn=m_Gmn[analysis_counter];
        nodesets=m_nodesets[analysis_counter];
        yg=m_yg[analysis_counter];
        ys=m_ys[analysis_counter];

        /*Now, plug analysis_counter and analysis_type inside the parameters: */
        this->parameters->SetParam(analysis_counter,AnalysisCounterEnum);
        this->parameters->SetParam(analysis_type,AnalysisTypeEnum);
        this->parameters->SetParam(configuration_type,ConfigurationTypeEnum);

        /*configure elements, loads and nodes, for this new analysis: */
        this->elements->  Configure(elements,loads, nodes,vertices, materials,parameters);
        this->nodes->     Configure(elements,loads, nodes,vertices, materials,parameters);
        this->loads->     Configure(elements, loads, nodes,vertices, materials,parameters);

}
/*}}}1*/
/*FUNCTION FemModel::SetCurrentConfiguration(int configuration_type){{{1*/
void FemModel::SetCurrentConfiguration(int configuration_type){

        /*overload: analysis_type = configuration_type: */
        this->SetCurrentConfiguration(configuration_type,configuration_type);
}
/*}}}1*/