Skip to content
Snippets Groups Projects
Commit 8633fa99 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

fix multistep dg output

parent 5113d0a2
No related branches found
No related tags found
No related merge requests found
...@@ -12,7 +12,7 @@ ...@@ -12,7 +12,7 @@
#include <sstream> #include <sstream>
#include "MElement.h" #include "MElement.h"
dgDofContainer::dgDofContainer (dgGroupCollection *groups, int nbFields): dgDofContainer::dgDofContainer (dgGroupCollection *groups, int nbFields):
_groups(*groups) _groups(*groups), _mshStep(0)
{ {
int _dataSize = 0; int _dataSize = 0;
_dataSizeGhost=0; _dataSizeGhost=0;
...@@ -213,11 +213,13 @@ void dgDofContainer::L2Projection(std::string functionName){ ...@@ -213,11 +213,13 @@ void dgDofContainer::L2Projection(std::string functionName){
} }
void dgDofContainer::exportMsh(const std::string name){ void dgDofContainer::exportMsh(const std::string name)
{
// the elementnodedata::export does not work !! // the elementnodedata::export does not work !!
for (int ICOMP = 0; ICOMP<_nbFields;++ICOMP){ for (int ICOMP = 0; ICOMP<_nbFields;++ICOMP){
std::ostringstream name_oss; std::ostringstream name_oss, name_view;
name_view<<"comp_"<<ICOMP;
name_oss<<name<<"_COMP_"<<ICOMP<<".msh"; name_oss<<name<<"_COMP_"<<ICOMP<<".msh";
if(Msg::GetCommSize()>1) if(Msg::GetCommSize()>1)
name_oss<<"_"<<Msg::GetCommRank(); name_oss<<"_"<<Msg::GetCommRank();
...@@ -229,11 +231,11 @@ void dgDofContainer::exportMsh(const std::string name){ ...@@ -229,11 +231,11 @@ void dgDofContainer::exportMsh(const std::string name){
fprintf(f,"$MeshFormat\n2.1 0 8\n$EndMeshFormat\n"); fprintf(f,"$MeshFormat\n2.1 0 8\n$EndMeshFormat\n");
fprintf(f,"$ElementNodeData\n"); fprintf(f,"$ElementNodeData\n");
fprintf(f,"1\n"); fprintf(f,"1\n");
fprintf(f,"\"%s\"\n",name.c_str()); fprintf(f,"\"%s\"\n",name_view.str().c_str());
fprintf(f,"1\n"); fprintf(f,"1\n");
fprintf(f,"0.0\n"); fprintf(f,"%d\n", _mshStep); // should print actual time here
fprintf(f,"%d\n", Msg::GetCommSize() > 1 ? 4 : 3); fprintf(f,"%d\n", Msg::GetCommSize() > 1 ? 4 : 3);
fprintf(f,"0\n 1\n %d\n",COUNT); fprintf(f,"%d\n 1\n %d\n", _mshStep, COUNT);
if(Msg::GetCommSize() > 1) fprintf(f,"%d\n", Msg::GetCommRank()); if(Msg::GetCommSize() > 1) fprintf(f,"%d\n", Msg::GetCommRank());
for (int i=0;i < _groups.getNbElementGroups() ;i++){ for (int i=0;i < _groups.getNbElementGroups() ;i++){
dgGroupOfElements *group = _groups.getElementGroup(i); dgGroupOfElements *group = _groups.getElementGroup(i);
...@@ -251,7 +253,9 @@ void dgDofContainer::exportMsh(const std::string name){ ...@@ -251,7 +253,9 @@ void dgDofContainer::exportMsh(const std::string name){
fprintf(f,"$EndElementNodeData\n"); fprintf(f,"$EndElementNodeData\n");
fclose(f); fclose(f);
} }
return;
_mshStep++;
// we should discuss that : we do a copy of the solution, this should // we should discuss that : we do a copy of the solution, this should
// be avoided ! // be avoided !
......
...@@ -23,6 +23,7 @@ private: ...@@ -23,6 +23,7 @@ private:
double *sendBuf, *recvBuf; double *sendBuf, *recvBuf;
std::vector<fullMatrix<double> *> _dataProxys; // proxys std::vector<fullMatrix<double> *> _dataProxys; // proxys
std::map<const dgGroupOfElements*,int> _groupId; std::map<const dgGroupOfElements*,int> _groupId;
int _mshStep;
public: public:
void scale(double f); void scale(double f);
double norm(); double norm();
......
...@@ -21,12 +21,10 @@ eigenSolver::eigenSolver(dofManager<double> *manager, std::string A, ...@@ -21,12 +21,10 @@ eigenSolver::eigenSolver(dofManager<double> *manager, std::string A,
_B = dynamic_cast<linearSystemPETSc<double>*>(manager->getLinearSystem(B)); _B = dynamic_cast<linearSystemPETSc<double>*>(manager->getLinearSystem(B));
if(!_B) Msg::Error("Could not find PETSc system '%s'", B.c_str()); if(!_B) Msg::Error("Could not find PETSc system '%s'", B.c_str());
} }
} }
bool eigenSolver::solve(int numEigenValues, std::string which) bool eigenSolver::solve(int numEigenValues, std::string which)
{ {
if(!_A) return false; if(!_A) return false;
Mat A = _A->getMatrix(); Mat A = _A->getMatrix();
Mat B = _B ? _B->getMatrix() : PETSC_NULL; Mat B = _B ? _B->getMatrix() : PETSC_NULL;
...@@ -111,11 +109,10 @@ bool eigenSolver::solve(int numEigenValues, std::string which) ...@@ -111,11 +109,10 @@ bool eigenSolver::solve(int numEigenValues, std::string which)
_try(EPSGetConverged(eps, &nconv)); _try(EPSGetConverged(eps, &nconv));
Msg::Info("SLEPc number of converged eigenpairs: %d", nconv); Msg::Info("SLEPc number of converged eigenpairs: %d", nconv);
if (nconv>0) {
// ignore additional eigenvalues if we get more than what we asked // ignore additional eigenvalues if we get more than what we asked
if(nconv > nev) nconv = nev; if(nconv > nev) nconv = nev;
if (nconv > 0) {
Vec xr, xi; Vec xr, xi;
_try(MatGetVecs(A, PETSC_NULL, &xr)); _try(MatGetVecs(A, PETSC_NULL, &xr));
_try(MatGetVecs(A, PETSC_NULL, &xi)); _try(MatGetVecs(A, PETSC_NULL, &xi));
...@@ -151,16 +148,12 @@ bool eigenSolver::solve(int numEigenValues, std::string which) ...@@ -151,16 +148,12 @@ bool eigenSolver::solve(int numEigenValues, std::string which)
} }
_eigenVectors.push_back(ev); _eigenVectors.push_back(ev);
} }
// cleanup
_try(EPSDestroy(eps));
_try(VecDestroy(xr)); _try(VecDestroy(xr));
_try(VecDestroy(xi)); _try(VecDestroy(xi));
_try(SlepcFinalize());
} }
_try(EPSDestroy(eps));
_try(SlepcFinalize());
if(reason == EPS_CONVERGED_TOL){ if(reason == EPS_CONVERGED_TOL){
Msg::Info("SLEPc done"); Msg::Info("SLEPc done");
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment