OilfieldToolsNet/IFEM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Added output of projected secondary solution variables to VTF

Browse Source 
git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@906 e10b68d5-8a6e-419e-a041-bce267b0401d
This commit is contained in:
kmo
2011-04-12 08:14:41 +00:00
committed by Knut Morten Okstad
parent 8e4aeca49f
commit 59d9fdc990
27 changed files with 213 additions and 100 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
src/ASM/ASMbase.C
View File

@@ -539,7 +539,7 @@ bool ASMbase::evalSolution (Matrix&, const Vector&,
}
bool ASMbase::evalSolution (Matrix&, const Integrand&, const int*) const
bool ASMbase::evalSolution (Matrix&, const Integrand&, const int*, bool) const
{
std::cerr <<" *** ASMBase::evalSolution: Must be implemented in sub-class."
<< std::endl;

5
src/ASM/ASMbase.h
View File

@@ -270,14 +270,17 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
//! \param[in] project Flag indicating if the projected solution is wanted
//!
//! \details The secondary solution is derived from the primary solution,
//! which is assumed to be stored within the \a integrand for current patch.
//! If \a npe is NULL, the solution is evaluated at the Greville points and
//! then projected onto the spline basis to obtain the control point values,
//! which then are returned through \a sField.
//! If \a npe is not NULL and \a project is \e true, the solution is also
//! projected onto the spline basis, and then evaluated at the \a npe points.
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe = 0) const;
const int* npe = 0, bool project = false) const;
//! \brief Evaluates the secondary solution field at the given points.
//! \param[out] sField Solution field

59
src/ASM/ASMs1D.C
View File

@@ -681,21 +681,24 @@ bool ASMs1D::tesselate (ElementBlock& grid, const int* npe) const
// Evaluate the spline curve at all points
size_t nx = gpar.size();
RealArray XYZ(3*nx,0.0);
// Establish the block grid coordinates
size_t i, j;
for (i = j = 0; i < nx; i++, j += 3)
RealArray XYZ(curv->dimension()*nx);
//curv->gridEvaluator(XYZ,gpar);
//TODO: Replace this loop by the above line when available in GoTools
size_t m, jp;
Go::Point pt;
for (m = jp = 0; m < nx; m++)
{
Go::Point pt;
curv->point(pt,gpar[i]);
curv->point(pt,gpar[m]);
for (int k = 0; k < pt.size(); k++)
XYZ[j+k] = pt[k];
XYZ[jp++] = pt[k];
}
// Establish the block grid coordinates
size_t i, j, l;
grid.resize(nx);
for (i = j = 0; i < grid.getNoNodes(); i++, j += 3)
grid.setCoor(i,XYZ[j],XYZ[j+1],XYZ[j+2]);
for (i = l = 0; i < grid.getNoNodes(); i++, l += curv->dimension())
for (j = 0; j < nsd; j++)
grid.setCoor(i,j,XYZ[l+j]);
// Establish the block grid topology
int n[2], ip = 0;
@@ -703,8 +706,8 @@ bool ASMs1D::tesselate (ElementBlock& grid, const int* npe) const
n[1] = n[0] + 1;
for (i = 1; i < nx; i++)
for (j = 0; j < 2; j++)
grid.setNode(ip++,n[j]++);
for (l = 0; l < 2; l++)
grid.setNode(ip++,n[l]++);
return true;
}
@@ -767,15 +770,38 @@ bool ASMs1D::evalSolution (Matrix& sField, const Vector& locSol,
bool ASMs1D::evalSolution (Matrix& sField, const Integrand& integrand,
const int* npe) const
const int* npe, bool project) const
{
if (npe)
{
// Compute parameter values of the result sampling points
RealArray gpar;
if (this->getGridParameters(gpar,npe[0]-1))
// Evaluate the secondary solution at all sampling points
return this->evalSolution(sField,integrand,&gpar);
if (project)
{
// Project the secondary solution onto the spline basis
Go::SplineCurve* c = this->projectSolution(integrand);
if (c)
{
// Evaluate the projected field at the result sampling points
//const Vector& svec = sField; // using utl::matrix cast operator
sField.resize(c->dimension(),gpar.size());
//c->gridEvaluator(svec,gpar);
//TODO: Replace this loop by the above line when available in GoTools
Go::Point pt;
for (size_t j = 0; j < sField.cols(); j++)
{
c->point(pt,gpar[j]);
for (size_t i = 0; i < sField.rows(); i++)
sField(i+1,j+1) = pt[i];
}
delete c;
return true;
}
}
else
// Evaluate the secondary solution directly at all sampling points
return this->evalSolution(sField,integrand,&gpar);
}
else
{
@@ -791,6 +817,7 @@ bool ASMs1D::evalSolution (Matrix& sField, const Integrand& integrand,
}
}
std::cerr <<" *** ASMs1D::evalSolution: Failure!"<< std::endl;
return false;
}
@@ -805,7 +832,7 @@ Go::SplineCurve* ASMs1D::projectSolution (const Integrand& integrand) const
{
//TODO: This requires the method Go::CurveInterpolator::regularInterpolation
// which is not yet present in GoTools.
std::cerr <<"ASMs1D::projectSolution: Not implemented yet!"<< std::endl;
std::cerr <<" *** ASMs1D::projectSolution: Not implemented yet!"<< std::endl;
return 0;
}

5
src/ASM/ASMs1D.h
View File

@@ -154,14 +154,17 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
//! \param[in] project Flag indicating if the projected solution is wanted
//!
//! \details The secondary solution is derived from the primary solution,
//! which is assumed to be stored within the \a integrand for current patch.
//! If \a npe is NULL, the solution is evaluated at the Greville points and
//! then projected onto the spline basis to obtain the control point values,
//! which then are returned through \a sField.
//! If \a npe is not NULL and \a project is \e true, the solution is also
//! projected onto the spline basis, and then evaluated at the \a npe points.
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe = 0) const;
const int* npe = 0, bool project = false) const;
//! \brief Projects the secondary solution field onto the primary basis.
//! \param[in] integrand Object with problem-specific data and methods

16
src/ASM/ASMs1DLag.C
View File

@@ -73,16 +73,16 @@ bool ASMs1DLag::generateFEMTopology ()
coord.resize(nx);
// Evaluate the nodal coordinates
Go::Point pt;
for (size_t i = 0; i < nx; i++)
{
Go::Point pt;
curv->point(pt,gpar[i]);
for (int k = 0; k < pt.size(); k++)
coord[i][k] = pt[k];
MLGN[i] = ++gNod;
}
// Connectivety array: local --> global node relation
// Connectivity array: local --> global node relation
MLGE.resize(nel);
MNPC.resize(nel);
@@ -299,8 +299,16 @@ int ASMs1DLag::evalPoint (const double* xi, double* param, Vec3& X) const
}
bool ASMs1DLag::tesselate (ElementBlock& grid, const int*) const
bool ASMs1DLag::tesselate (ElementBlock& grid, const int* npe) const
{
const int p1 = curv->order();
if (p1 != npe[0])
{
std::cout <<"\nLagrange elements: The number of visualization points is "
<< p1 <<" by default\n"<< std::endl;
const_cast<int*>(npe)[0] = p1;
}
size_t i, l;
grid.resize(nx);
@@ -348,7 +356,7 @@ bool ASMs1DLag::evalSolution (Matrix&, const Vector&,
bool ASMs1DLag::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

9
src/ASM/ASMs1DLag.h
View File

@@ -39,7 +39,7 @@ public:
//! \brief Generates the finite element topology data for the patch.
//! \details The data generated are the element-to-node connectivity array,
//! and the global node and element numbers.
//! the nodal coordinate array, as well as global node and element numbers.
virtual bool generateFEMTopology();
//! \brief Clears the contents of the patch, making it empty.
@@ -86,8 +86,9 @@ public:
//! \brief Creates a line element model of this patch for visualization.
//! \param[out] grid The generated line grid
//! \param[in] npe Number of visualization nodes over each knot span
//! \note The number of element nodes must be set in \a grid on input.
virtual bool tesselate(ElementBlock& grid, const int*) const;
virtual bool tesselate(ElementBlock& grid, const int* npe) const;
//! \brief Evaluates the primary solution field at all visualization points.
//! \details The number of visualization points is the same as the order of
@@ -110,7 +111,7 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int*) const;
const int*, bool = false) const;
//! \brief Evaluates the secondary solution field at the given points.
//! \param[out] sField Solution field
@@ -135,7 +136,7 @@ protected:
virtual void getNodalCoordinates(Matrix& X) const;
//! \brief Returns the number of nodal points in the patch.
virtual int getSize() const { return nx; }
virtual int getSize(int = 0) const { return nx; }
private:
//! \brief The number of nodes in each direction for the patch

2
src/ASM/ASMs1DSpec.C
View File

@@ -157,7 +157,7 @@ bool ASMs1DSpec::integrate (Integrand& integrand, int lIndex,
bool ASMs1DSpec::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);
if (!curv) return false;

3
src/ASM/ASMs1DSpec.h
View File

@@ -65,9 +65,8 @@ public:
//! \brief Evaluates the secondary solution field at all visualization points.
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe) const;
const int*, bool = false) const;
protected:

45
src/ASM/ASMs2D.C
View File

@@ -1178,14 +1178,10 @@ bool ASMs2D::tesselate (ElementBlock& grid, const int* npe) const
// Establish the block grid coordinates
size_t i, j, l;
double X[3] = { 0.0, 0.0, 0.0 };
grid.resize(nx,ny);
for (i = j = 0; i < grid.getNoNodes(); i++, j += surf->dimension())
{
for (l = 0; l < nsd; l++)
X[l] = XYZ[j+l];
grid.setCoor(i,X[0],X[1],X[2]);
}
for (i = l = 0; i < grid.getNoNodes(); i++, l += surf->dimension())
for (j = 0; j < nsd; j++)
grid.setCoor(i,j,XYZ[l+j]);
// Establish the block grid topology
int n[4], ip = 0;
@@ -1288,20 +1284,31 @@ bool ASMs2D::evalSolution (Matrix& sField, const Vector& locSol,
bool ASMs2D::evalSolution (Matrix& sField, const Integrand& integrand,
const int* npe) const
const int* npe, bool project) const
{
bool retVal = true;
if (npe)
{
// Compute parameter values of the result sampling points
RealArray gpar[2];
for (int dir = 0; dir < 2 && retVal; dir++)
retVal = this->getGridParameters(gpar[dir],dir,npe[dir]-1);
// Evaluate the secondary solution at all sampling points
if (retVal)
retVal = this->evalSolution(sField,integrand,gpar);
if (this->getGridParameters(gpar[0],0,npe[0]-1) &&
this->getGridParameters(gpar[1],1,npe[1]-1))
if (project)
{
// Project the secondary solution onto the spline basis
Go::SplineSurface* s = this->projectSolution(integrand);
if (s)
{
// Evaluate the projected field at the result sampling points
const Vector& svec = sField; // using utl::matrix cast operator
sField.resize(s->dimension(),gpar[0].size()*gpar[1].size());
s->gridEvaluator(const_cast<Vector&>(svec),gpar[0],gpar[1]);
delete s;
return true;
}
}
else
// Evaluate the secondary solution directly at all sampling points
return this->evalSolution(sField,integrand,gpar);
}
else
{
@@ -1313,12 +1320,12 @@ bool ASMs2D::evalSolution (Matrix& sField, const Integrand& integrand,
sField.resize(s->dimension(),s->numCoefs_u()*s->numCoefs_v());
sField.fill(&(*s->coefs_begin()));
delete s;
return true;
}
else
retVal = false;
}
return retVal;
std::cerr <<" *** ASMs2D::evalSolution: Failure!"<< std::endl;
return false;
}

5
src/ASM/ASMs2D.h
View File

@@ -235,14 +235,17 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
//! \param[in] project Flag indicating if the projected solution is wanted
//!
//! \details The secondary solution is derived from the primary solution,
//! which is assumed to be stored within the \a integrand for current patch.
//! If \a npe is NULL, the solution is evaluated at the Greville points and
//! then projected onto the spline basis to obtain the control point values,
//! which then are returned through \a sField.
//! If \a npe is not NULL and \a project is \e true, the solution is also
//! projected onto the spline basis, and then evaluated at the \a npe points
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe = 0) const;
const int* npe = 0, bool project = false) const;
//! \brief Projects the secondary solution field onto the primary basis.
//! \param[in] integrand Object with problem-specific data and methods

4
src/ASM/ASMs2DLag.C
View File

@@ -93,7 +93,7 @@ bool ASMs2DLag::generateFEMTopology ()
// Number of nodes per element
const int nen = p1*p2;
// Connectivety array: local --> global node relation
// Connectivity array: local --> global node relation
MLGE.resize(nel);
MNPC.resize(nel);
@@ -453,7 +453,7 @@ bool ASMs2DLag::evalSolution (Matrix&, const Vector&,
bool ASMs2DLag::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

2
src/ASM/ASMs2DLag.h
View File

@@ -112,7 +112,7 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int*) const;
const int*, bool = false) const;
//! \brief Evaluates the secondary solution field at the given points.
//! \param[out] sField Solution field

2
src/ASM/ASMs2DSpec.C
View File

@@ -264,7 +264,7 @@ bool ASMs2DSpec::integrate (Integrand& integrand, int lIndex,
bool ASMs2DSpec::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

3
src/ASM/ASMs2DSpec.h
View File

@@ -65,9 +65,8 @@ public:
//! \brief Evaluates the secondary solution field at all visualization points.
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe) const;
const int*, bool = false) const;
protected:

4
src/ASM/ASMs2DmxLag.C
View File

@@ -126,7 +126,7 @@ bool ASMs2DmxLag::generateFEMTopology ()
const int nelx = (nx2-1)/(p1-1);
const int nely = (ny2-1)/(p2-1);
// Add connectivety for second basis: local --> global node relation
// Add connectivity for second basis: local --> global node relation
int i, j, a, b, iel;
for (j = iel = 0; j < nely; j++)
for (i = 0; i < nelx; i++, iel++)
@@ -414,7 +414,7 @@ bool ASMs2DmxLag::evalSolution (Matrix& sField, const Vector& locSol,
bool ASMs2DmxLag::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

2
src/ASM/ASMs2DmxLag.h
View File

@@ -118,7 +118,7 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int*) const;
const int*, bool = false) const;
//! \brief Extracts nodal results for this patch from the global vector.
//! \param[in] globVec Global solution vector in DOF-order

37
src/ASM/ASMs3D.C
View File

@@ -1821,20 +1821,33 @@ bool ASMs3D::evalSolution (Matrix& sField, const Vector& locSol,
bool ASMs3D::evalSolution (Matrix& sField, const Integrand& integrand,
const int* npe) const
const int* npe, bool project) const
{
bool retVal = true;
if (npe)
{
// Compute parameter values of the result sampling points
RealArray gpar[3];
for (int dir = 0; dir < 3 && retVal; dir++)
retVal = this->getGridParameters(gpar[dir],dir,npe[dir]-1);
// Evaluate the secondary solution at all sampling points
if (retVal)
retVal = this->evalSolution(sField,integrand,gpar);
if (this->getGridParameters(gpar[0],0,npe[0]-1) &&
this->getGridParameters(gpar[1],1,npe[1]-1) &&
this->getGridParameters(gpar[2],2,npe[2]-1))
if (project)
{
// Project the secondary solution onto the spline basis
Go::SplineVolume* v = this->projectSolution(integrand);
if (v)
{
// Evaluate the projected field at the result sampling points
const Vector& svec = sField; // using utl::matrix cast operator
sField.resize(v->dimension(),
gpar[0].size()*gpar[1].size()*gpar[2].size());
v->gridEvaluator(const_cast<Vector&>(svec),gpar[0],gpar[1],gpar[2]);
delete v;
return true;
}
}
else
// Evaluate the secondary solution at all sampling points
return this->evalSolution(sField,integrand,gpar);
}
else
{
@@ -1847,12 +1860,12 @@ bool ASMs3D::evalSolution (Matrix& sField, const Integrand& integrand,
v->numCoefs(0)*v->numCoefs(1)*v->numCoefs(2));
sField.fill(&(*v->coefs_begin()));
delete v;
return true;
}
else
retVal = false;
}
return retVal;
std::cerr <<" *** ASMs3D::evalSolution: Failure!"<< std::endl;
return false;
}

5
src/ASM/ASMs3D.h
View File

@@ -300,14 +300,17 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
//! \param[in] project Flag indicating if the projected solution is wanted
//!
//! \details The secondary solution is derived from the primary solution,
//! which is assumed to be stored within the \a integrand for current patch.
//! If \a npe is NULL, the solution is evaluated at the Greville points and
//! then projected onto the spline basis to obtain the control point values,
//! which then are returned through \a sField.
//! If \a npe is not NULL and \a project is \e true, the solution is also
//! projected onto the spline basis, and then evaluated at the \a npe points
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe = 0) const;
const int* npe = 0, bool project = false) const;
//! \brief Projects the secondary solution field onto the primary basis.
//! \param[in] integrand Object with problem-specific data and methods

4
src/ASM/ASMs3DLag.C
View File

@@ -98,7 +98,7 @@ bool ASMs3DLag::generateFEMTopology()
// Number of nodes in a xy-surface of an element
const int ct = p1*p2;
// Connectivety array: local --> global node relation
// Connectivity array: local --> global node relation
MLGE.resize(nel);
MNPC.resize(nel);
@@ -614,7 +614,7 @@ bool ASMs3DLag::evalSolution (Matrix&, const Vector&,
bool ASMs3DLag::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

2
src/ASM/ASMs3DLag.h
View File

@@ -120,7 +120,7 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int*) const;
const int*, bool = false) const;
//! \brief Evaluates the secondary solution field at the given points.
//! \param[out] sField Solution field

2
src/ASM/ASMs3DSpec.C
View File

@@ -408,7 +408,7 @@ bool ASMs3DSpec::integrateEdge (Integrand& integrand, int lEdge,
bool ASMs3DSpec::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

3
src/ASM/ASMs3DSpec.h
View File

@@ -73,9 +73,8 @@ public:
//! \brief Evaluates the secondary solution field at all visualization points.
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
//! \param[in] npe Number of visualization nodes over each knot span
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int* npe) const;
const int*, bool = false) const;
protected:

4
src/ASM/ASMs3DmxLag.C
View File

@@ -131,7 +131,7 @@ bool ASMs3DmxLag::generateFEMTopology ()
const int nely = (ny2-1)/(p2-1);
const int nelz = (nz2-1)/(p3-1);
// Add connectivety for second basis: local --> global node relation
// Add connectivity for second basis: local --> global node relation
int i, j, k, a, b, c, iel;
for (k = iel = 0; k < nelz; k++)
for (j = 0; j < nely; j++)
@@ -451,7 +451,7 @@ bool ASMs3DmxLag::evalSolution (Matrix& sField, const Vector& locSol,
bool ASMs3DmxLag::evalSolution (Matrix& sField, const Integrand& integrand,
const int*) const
const int*, bool) const
{
sField.resize(0,0);

2
src/ASM/ASMs3DmxLag.h
View File

@@ -118,7 +118,7 @@ public:
//! \param[out] sField Solution field
//! \param[in] integrand Object with problem-specific data and methods
virtual bool evalSolution(Matrix& sField, const Integrand& integrand,
const int*) const;
const int*, bool = false) const;
//! \brief Extracts nodal results for this patch from the global vector.
//! \param[in] globVec Global solution vector in DOF-order

71
src/SIM/SIMbase.C
View File

@@ -926,6 +926,7 @@ bool SIMbase::writeGlvBC (const int* nViz, int& nBlock) const
if (msgLevel > 1)
std::cout <<"Writing boundary conditions for patch "<< i+1 << std::endl;
if (!myModel[i]->evalSolution(field,bc,nViz))
return false;
@@ -980,7 +981,7 @@ bool SIMbase::writeGlvV (const Vector& vec, const char* fieldName,
Matrix field;
Vector lovec;
int geomID = 0, idBlock = 2;
int geomID = 0;
std::vector<int> vID;
for (size_t i = 0; i < myModel.size(); i++)
@@ -989,6 +990,7 @@ bool SIMbase::writeGlvV (const Vector& vec, const char* fieldName,
if (msgLevel > 1)
std::cout <<"Writing vector field for patch "<< i+1 << std::endl;
myModel[i]->extractNodeVec(vec,lovec);
if (!myModel[i]->evalSolution(field,lovec,nViz))
return false;
@@ -999,6 +1001,7 @@ bool SIMbase::writeGlvV (const Vector& vec, const char* fieldName,
vID.push_back(nBlock);
}
int idBlock = 2;
return myVtf->writeVblk(vID,fieldName,idBlock,iStep);
}
@@ -1016,8 +1019,9 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
Matrix field;
size_t i, j;
int geomID = 0, idBlock = 10;
std::vector<int> vID, dID[14], sID[14];
int geomID = 0;
const size_t sMAX = 21;
std::vector<int> vID, dID[3], sID[sMAX];
bool haveAsol = false;
bool scalarEq = myModel.front()->getNoFields() == 1;
if (mySol)
@@ -1032,6 +1036,9 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
if (msgLevel > 1)
std::cout <<"Writing FE solution for patch "<< i+1 << std::endl;
// 1. Evaluate primary solution variables
myModel[i]->extractNodeVec(psol,myProblem->getSolution());
if (!myModel[i]->evalSolution(field,myProblem->getSolution(),nViz))
return false;
@@ -1059,6 +1066,8 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
if (psolOnly) continue; // skip secondary solution
// 2. Direct evaluation of secondary solution variables
LocalSystem::patch = i;
if (!myModel[i]->evalSolution(field,*myProblem,nViz))
return false;
@@ -1070,16 +1079,33 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
vID.push_back(nBlock);
size_t k = 0;
for (j = 1; j <= field.rows() && k < 14; j++)
for (j = 1; j <= field.rows() && k < sMAX; j++)
if (!myVtf->writeNres(field.getRow(j),++nBlock,geomID))
return false;
else
sID[k++].push_back(nBlock);
if (discretization == Spline)
{
// 3. Projection of secondary solution variables (splines only)
if (!myModel[i]->evalSolution(field,*myProblem,nViz,true))
return false;
for (j = 1; j <= field.rows() && k < sMAX; j++)
if (!myVtf->writeNres(field.getRow(j),++nBlock,geomID))
return false;
else
sID[k++].push_back(nBlock);
}
if (haveAsol)
{
// 4. Evaluate analytical solution variables
if (msgLevel > 1)
std::cout <<"Writing exact solution for patch "<< i+1 << std::endl;
const ElementBlock* grid = myVtf->getBlock(geomID);
std::vector<Vec3>::const_iterator cit = grid->begin_XYZ();
Vector solPt; field.fill(0.0);
@@ -1095,7 +1121,8 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
if (ok)
field.fillColumn(j,solPt);
}
for (j = 1; j <= field.rows() && k < 14 && ok; j++)
for (j = 1; j <= field.rows() && k < sMAX && ok; j++)
if (!myVtf->writeNres(field.getRow(j),++nBlock,geomID))
return false;
else
@@ -1103,15 +1130,17 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
}
}
// Write result block identifications
int idBlock = 10;
if (scalarEq)
{
if (!psolOnly)
if (!myVtf->writeVblk(vID,"Flux",idBlock,iStep))
return false;
if (!myVtf->writeSblk(dID[0],"u",++idBlock,iStep))
return false;
if (psolOnly) return true;
if (!myVtf->writeVblk(vID,"Flux",idBlock,iStep))
return false;
}
else
{
@@ -1129,13 +1158,18 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
if (psolOnly) return true;
size_t nf = myProblem->getNoFields();
for (j = 0; j < nf && !sID[j].empty(); j++)
if (!myVtf->writeSblk(sID[j],myProblem->getFieldLabel(j,haveAsol?"FE":0),
for (i = j = 0; i < nf && !sID[j].empty(); i++, j++)
if (!myVtf->writeSblk(sID[j],myProblem->getFieldLabel(i,haveAsol?"FE":0),
++idBlock,iStep)) return false;
if (discretization == Spline)
for (i = 0; i < nf && !sID[j].empty(); i++, j++)
if (!myVtf->writeSblk(sID[j],myProblem->getFieldLabel(i,"Projected"),
++idBlock,iStep)) return false;
if (haveAsol)
for (j = 0; j < nf && !sID[nf+j].empty(); j++)
if (!myVtf->writeSblk(sID[nf+j],myProblem->getFieldLabel(j,"Exact"),
for (i = 0; i < nf && !sID[j].empty(); i++, j++)
if (!myVtf->writeSblk(sID[j],myProblem->getFieldLabel(i,"Exact"),
++idBlock,iStep)) return false;
return true;
@@ -1164,7 +1198,7 @@ bool SIMbase::writeGlvM (const Mode& mode, bool freq,
Vector displ;
Matrix field;
int geomID = 0, idBlock = 10;
int geomID = 0;
std::vector<int> vID;
for (size_t i = 0; i < myModel.size(); i++)
@@ -1184,6 +1218,7 @@ bool SIMbase::writeGlvM (const Mode& mode, bool freq,
}
if (myModel.size() > 1 && msgLevel > 1) std::cout << std::endl;
int idBlock = 10;
if (!myVtf->writeDblk(vID,"Mode Shape",idBlock,mode.eigNo))
return false;
@@ -1200,7 +1235,7 @@ bool SIMbase::writeGlvN (const Matrix& norms, int iStep, int& nBlock)
return false;
Matrix field;
int idBlock = 100, geomID = 0;
int geomID = 0;
std::vector<int> sID[3];
size_t i, j, k;
@@ -1208,9 +1243,10 @@ bool SIMbase::writeGlvN (const Matrix& norms, int iStep, int& nBlock)
{
if (myModel[i]->empty()) continue; // skip empty patches
geomID++;
if (msgLevel > 1)
std::cout <<"Writing element norms for patch "<< i+1 << std::endl;
geomID++;
myModel[i]->extractElmRes(norms,field);
for (j = k = 0; j < field.rows() && j < 4; j++)
@@ -1227,6 +1263,7 @@ bool SIMbase::writeGlvN (const Matrix& norms, int iStep, int& nBlock)
"a(e,e)^0.5, e=u-u^h"
};
int idBlock = 100;
for (j = 0; j < 3; j++)
if (!sID[j].empty())
if (!myVtf->writeSblk(sID[j],label[j],++idBlock,iStep))

11
src/Utility/ElementBlock.C
View File

@@ -1,4 +1,4 @@
// $Id: ElementBlock.C,v 1.4 2010-09-05 11:29:18 kmo Exp $
// $Id$
//==============================================================================
//!
//! \file ElementBlock.C
@@ -53,6 +53,15 @@ bool ElementBlock::setCoor (size_t i, real x, real y, real z)
}
bool ElementBlock::setCoor (size_t i, size_t j, real x)
{
if (i >= coord.size() || j >= 3) return false;
coord[i][j] = x;
return true;
}
bool ElementBlock::setNode (size_t i, int nodeNumb)
{
if (i >= MMNPC.size()) return false;

4
src/Utility/ElementBlock.h
View File

@@ -1,4 +1,4 @@
// $Id: ElementBlock.h,v 1.2 2010-03-10 13:59:02 kmo Exp $
// $Id$
//==============================================================================
//!
//! \file ElementBlock.h
@@ -35,6 +35,8 @@ public:
//! \brief Defines the coordinates of node \a i
bool setCoor(size_t i, real x, real y, real z);
//! \brief Defines the \a j'th coordinate of node \a i
bool setCoor(size_t i, size_t j, real x);
//! \brief Defines the global number of element node \a i
bool setNode(size_t i, int nodeNumb);