From 29ac626f987a10c39693672bbbc6f462845a178e Mon Sep 17 00:00:00 2001
From: Gauthier Becker <gauthierbecker@gmail.com>
Date: Fri, 25 Mar 2011 15:26:45 +0000
Subject: [PATCH] Fix bug in NLS: creation of interface ans interdomain, Now no
Dof type suposition in solver Update all benchmark lua --> python
---
.../BoundaryConditions/nonLinearBC.h | 4 +-
NonLinearSolver/Domain/partDomain.h | 4 -
NonLinearSolver/Interface/MInterfaceElement.h | 2 +-
NonLinearSolver/Interface/MInterfaceLine.cpp | 3 +-
NonLinearSolver/Interface/MInterfaceLine.h | 5 +-
NonLinearSolver/contact/contactDomain.cpp | 18 +--
NonLinearSolver/contact/contactDomain.h | 24 ++-
.../contact/contactFunctionSpace.h | 44 ++----
NonLinearSolver/contact/contactTerms.h | 4 +-
.../contact/rigidCylinderContactTerms.cpp | 8 +-
.../contact/rigidCylinderContactTerms.h | 18 +--
NonLinearSolver/internalPoints/ipField.h | 17 +--
NonLinearSolver/internalPoints/ipstate.cpp | 4 +-
NonLinearSolver/nlsolver/CMakeLists.txt | 1 -
NonLinearSolver/nlsolver/Dof3IntType.h | 32 ----
NonLinearSolver/nlsolver/nlsolAlgorithms.h | 6 +-
.../nlsolver/nonLinearMechSolver.cpp | 26 ++--
.../nlsolver/nonLinearMechSolver.h | 11 +-
.../nlsolver/unknownDynamicField.cpp | 8 +-
.../nlsolver/unknownDynamicField.h | 2 +-
NonLinearSolver/nlsolver/unknownField.cpp | 98 ++++++++++--
NonLinearSolver/nlsolver/unknownField.h | 10 +-
.../nlsolver/unknownStaticField.cpp | 118 +-------------
NonLinearSolver/nlsolver/unknownStaticField.h | 2 +-
Post/PView.h | 18 +--
Post/PViewVertexArrays.cpp | 144 +++++++++---------
26 files changed, 261 insertions(+), 370 deletions(-)
delete mode 100644 NonLinearSolver/nlsolver/Dof3IntType.h
diff --git a/NonLinearSolver/BoundaryConditions/nonLinearBC.h b/NonLinearSolver/BoundaryConditions/nonLinearBC.h
index 835cdd359f..0d9873b8f6 100644
--- a/NonLinearSolver/BoundaryConditions/nonLinearBC.h
+++ b/NonLinearSolver/BoundaryConditions/nonLinearBC.h
@@ -56,7 +56,7 @@ class nonLinearDirichletBC : public nonLinearBoundaryCondition
class rigidContactBC : public nonLinearBoundaryCondition
{
public:
- rigidContactSpace *space;
+ rigidContactSpaceBase *space;
int _comp; // component
simpleFunctionTime<double> _f;
rigidContactBC(const int physMaster) : nonLinearBoundaryCondition(){
@@ -70,7 +70,7 @@ class rigidContactBC : public nonLinearBoundaryCondition
_f = source._f;
}
~rigidContactBC(){}
- void setSpace(rigidContactSpace *sp){space = sp;}
+ void setSpace(rigidContactSpaceBase *sp){space = sp;}
};
#endif
class nonLinearNeumannBC : public nonLinearBoundaryCondition
diff --git a/NonLinearSolver/Domain/partDomain.h b/NonLinearSolver/Domain/partDomain.h
index a703dd02de..57d5a2c80c 100644
--- a/NonLinearSolver/Domain/partDomain.h
+++ b/NonLinearSolver/Domain/partDomain.h
@@ -116,10 +116,6 @@ class partDomain{
// static void registerBindings(binding *b);
virtual void computeIPVariable(AllIPState *aips,const unknownField *ufield,const IPStateBase::whichState ws)=0;
- virtual void computeIpv(AllIPState *aips,MInterfaceElement *ie, IntPt *GP,const IPStateBase::whichState ws,
- partDomain* efMinus, partDomain *efPlus,materialLaw *mlawminus,
- materialLaw *mlawplus,fullVector<double> &dispm,
- fullVector<double> &dispp,const bool virt,const bool checkfrac=true)=0;
virtual void computeIpv(AllIPState *aips,MElement *e, IPStateBase::whichState ws,
materialLaw *mlaw,fullVector<double> &disp)=0;
virtual void setGaussIntegrationRule()=0;
diff --git a/NonLinearSolver/Interface/MInterfaceElement.h b/NonLinearSolver/Interface/MInterfaceElement.h
index dee64bd16a..1bc3d903c8 100644
--- a/NonLinearSolver/Interface/MInterfaceElement.h
+++ b/NonLinearSolver/Interface/MInterfaceElement.h
@@ -24,7 +24,7 @@ class MInterfaceElement{
virtual int getEdgeOrFaceNumber(const int index) const=0;
virtual bool isSameDirection(const int index) const=0;
// should return the element number !!
- virtual int getNum() const=0;
+ virtual int getNumber() const=0; // {return{this->getNum();} in your derived class it derived from an MElement* too !!
// compute the characteritic size of one element (This function can be defined as a method of MElement) ??
static double characSize(MElement *e); // Area/perimeter
};
diff --git a/NonLinearSolver/Interface/MInterfaceLine.cpp b/NonLinearSolver/Interface/MInterfaceLine.cpp
index a9eeed4903..d015455fc4 100644
--- a/NonLinearSolver/Interface/MInterfaceLine.cpp
+++ b/NonLinearSolver/Interface/MInterfaceLine.cpp
@@ -13,8 +13,7 @@
#include "MInterfaceLine.h"
MInterfaceLine::MInterfaceLine(std::vector<MVertex*> &v, int num, int part,
- MElement *e_minus, MElement *e_plus) : MLineN(v, num, part), MInterfaceElement(),
- _num(MElement::getGlobalNumber()) //avoid this ??
+ MElement *e_minus, MElement *e_plus) : MLineN(v, num, part), MInterfaceElement()
{
_numElem[0]=e_minus;
_numElem[1]=e_plus;
diff --git a/NonLinearSolver/Interface/MInterfaceLine.h b/NonLinearSolver/Interface/MInterfaceLine.h
index cd97f5b31d..46b91b7bc1 100644
--- a/NonLinearSolver/Interface/MInterfaceLine.h
+++ b/NonLinearSolver/Interface/MInterfaceLine.h
@@ -18,8 +18,6 @@
#include "MInterfaceElement.h"
class MInterfaceLine : public MLineN, public MInterfaceElement{ // or don't derivate but in this case a vector with the vertices of interface element has to be save ??
protected :
- // NUMBER DUPLICATION HOW TO AVOID THIS (MInterfaceElement : public MELelement but 2 MElement at this time and MELement has to be a virtual class ??)
- int _num;
// table of pointer on the two elements linked to the interface element
MElement *_numElem[2];
// edge's number linked to interface element of minus and plus element
@@ -35,8 +33,7 @@ class MInterfaceLine : public MLineN, public MInterfaceElement{ // or don't deri
~MInterfaceLine(){}
// Try to avoid this HOW
-// int getNum() const{MLineN::getNum();}
- int getNum() const{return _num;}
+ int getNumber() const{return this->getNum();}
// Give the number of minus 0 or plus 1 element
MElement* getElem(int index) const {return _numElem[index];}
diff --git a/NonLinearSolver/contact/contactDomain.cpp b/NonLinearSolver/contact/contactDomain.cpp
index 924d741534..044ffc4ad3 100644
--- a/NonLinearSolver/contact/contactDomain.cpp
+++ b/NonLinearSolver/contact/contactDomain.cpp
@@ -47,11 +47,11 @@ MElement* contactDomain::getFirstElement() const {
return (*it);
}
-rigidCylinderContact::rigidCylinderContact(const int tag, const int physMaster, const int physSlave, const int physPt1,
+rigidCylinderContactDomain::rigidCylinderContactDomain(const int tag, const int physMaster, const int physSlave, const int physPt1,
const int physPt2, const double penalty,
- const double h) : contactDomain(tag,physMaster,physSlave,
+ const double h,const double rho) : contactDomain(tag,physMaster,physSlave,
penalty,rigidCylinder,true),
- _thick(h){
+ _thick(h), _density(rho){
// void gauss integration
_integ = new QuadratureVoid();
@@ -98,18 +98,10 @@ rigidCylinderContact::rigidCylinderContact(const int tag, const int physMaster,
_axisDirector->normalize();
}
-void rigidCylinderContact::initializeTerms(const unknownField *ufield){
- rigidContactSpace *sp = dynamic_cast<rigidContactSpace*>(_space);
+void rigidCylinderContactDomain::initializeTerms(const unknownField *ufield){
+ rigidContactSpaceBase *sp = dynamic_cast<rigidContactSpaceBase*>(_space);
_massterm = new massRigidCylinder(this,sp);
_lterm = new forceRigidCylinderContact(this,sp,_thickContact,ufield);
rigidContactLinearTermBase<double> *rlterm = dynamic_cast<rigidContactLinearTermBase<double>*>(_lterm);
_bterm = new stiffnessRigidCylinderContact(sp,rlterm,ufield);
}
-
-void rigidCylinderContact::setDomain(partDomain *dom){
- _dom = dom;
- FunctionSpace<double> *spdom = dynamic_cast<FunctionSpace<double>*>(dom->getFunctionSpace());
- if(_space != NULL) delete _space;
- _space = new rigidContactSpace(_phys,spdom,_vergc);
-}
-
diff --git a/NonLinearSolver/contact/contactDomain.h b/NonLinearSolver/contact/contactDomain.h
index 885687a419..b48d0282d3 100644
--- a/NonLinearSolver/contact/contactDomain.h
+++ b/NonLinearSolver/contact/contactDomain.h
@@ -52,7 +52,6 @@ class contactDomain{
virtual void setTag(const int t){ _tag =t;}
virtual void setPhys(const int p){_phys =p;}
virtual void setPenalty(const double pe){_penalty=pe;}
- virtual void setDomain(partDomain *dom)=0;
virtual void setContactType(const int ct);
virtual void setContactType(const contact ct){_contype =ct;}
virtual int getTag() const{return _tag;}
@@ -69,11 +68,13 @@ class contactDomain{
virtual FunctionSpaceBase* getSpace() {return _space;}
virtual const FunctionSpaceBase* getSpace() const {return _space;}
virtual QuadratureBase* getGaussIntegration() const{return _integ;}
+ virtual void setDomainAndFunctionSpace(partDomain *dom)=0;
+ virtual void initializeTerms(const unknownField *ufield)=0;
// static void registerBindings(binding *b);
#endif
};
-class rigidCylinderContact : public contactDomain{
+class rigidCylinderContactDomain : public contactDomain{
protected:
double _length; // length of cylinder;
double _radius; // outer radius of cylinder;
@@ -82,27 +83,24 @@ class rigidCylinderContact : public contactDomain{
double _density; // density of cylinder Not a material law for now
MVertex *_vergc; // vertex of gravity center
SVector3 *_axisDirector; // normalized vector director of cylinder's axis
+#ifndef SWIG
public:
- rigidCylinderContact() : contactDomain(), _length(0.), _radius(0.), _thick(0.){
+ rigidCylinderContactDomain() : contactDomain(), _length(0.), _radius(0.), _thick(0.){
_contype = rigidCylinder;
_rigid = true;
_integ = new QuadratureVoid();
}
- rigidCylinderContact(const int tag, const int physMaster, const int physSlave, const int physpt1,
- const int physpt2,const double penalty,const double h);
- ~rigidCylinderContact(){delete _axisDirector;}
- void density(const double rho){ _density = rho;}
-// void setThicknessContact(const double h){_thickContact = h;}
-#ifndef SWIG
+ rigidCylinderContactDomain(const int tag, const int physMaster, const int physSlave, const int physpt1,
+ const int physpt2,const double penalty,const double h,const double rho);
+ ~rigidCylinderContactDomain(){delete _axisDirector;}
+ virtual void initializeTerms(const unknownField *ufield);
+ SVector3* getAxisDirector() const{return _axisDirector;}
+ virtual void setDomainAndFunctionSpace(partDomain *dom)=0;
virtual MVertex* getGC() const{return _vergc;}
double getDensity() const{return _density;}
double getLength() const{return _length;}
double getRadius() const{return _radius;}
double getThickness() const{return _thick;}
- void initializeTerms(const unknownField *ufield);
- SVector3* getAxisDirector() const{return _axisDirector;}
- virtual void setDomain(partDomain *dom);
#endif
-
};
#endif // CONTACTDOMAIN_H_
diff --git a/NonLinearSolver/contact/contactFunctionSpace.h b/NonLinearSolver/contact/contactFunctionSpace.h
index e2b67f9df9..ac510af28f 100644
--- a/NonLinearSolver/contact/contactFunctionSpace.h
+++ b/NonLinearSolver/contact/contactFunctionSpace.h
@@ -12,7 +12,6 @@
#ifndef CONTACTSPACE_H_
#define CONTACTSPACE_H_
#include "functionSpace.h"
-#include "Dof3IntType.h"
template <class T1> class ContactSpaceBase : public FunctionSpaceBase{
protected:
FunctionSpace<T1> *_space; // functional space of domain where contact is computed
@@ -42,44 +41,19 @@ template <class T1> class ContactSpaceBase : public FunctionSpaceBase{
}
};
-class rigidContactSpace : public ContactSpaceBase<double>{
+class rigidContactSpaceBase : public ContactSpaceBase<double>{
protected:
const MVertex *_gc; // Node of gravity center
public:
- rigidContactSpace(const int id, FunctionSpace<double> *sp, MVertex *ver) : ContactSpaceBase<double>(id,sp), _gc(ver){}
- rigidContactSpace(const int id, FunctionSpace<double> *sp, MVertex *ver, const int comp1,
+ rigidContactSpaceBase(const int id, FunctionSpace<double> *sp, MVertex *ver) : ContactSpaceBase<double>(id,sp), _gc(ver){}
+ rigidContactSpaceBase(const int id, FunctionSpace<double> *sp, MVertex *ver, const int comp1,
const int comp2 = -1, const int comp3 =-1) : ContactSpaceBase<double>(id,sp,comp1,comp2,comp3){}
- virtual void getKeys(MElement *ele,std::vector<Dof> &keys){
- _space->getKeys(ele,keys);
- for(int i=0;i<_comp.size();i++){
- keys.push_back(Dof(_gc->getNum(),Dof3IntType::createTypeWithThreeInts(_comp[i],_id)));
- }
- }
- virtual int getNumKeys(MElement *ele){
- int nkeysele = _space->getNumKeys(ele);
- return nkeysele + _comp.size();
- }
- virtual void getKeysOfGravityCenter(std::vector<Dof> &keys){
- for(int i=0;i<_comp.size(); i++)
- keys.push_back(Dof(_gc->getNum(),Dof3IntType::createTypeWithThreeInts(_comp[i],_id)));
- }
- virtual int getNumKeysOfGravityCenter(){
- return _comp.size();
- }
- virtual void getKeys(MElement *ele, const int ind, std::vector<Dof> &keys){
- // generate keys of element and select the good ones after LAGRANGE OK ?? CHANGE THIS HOW TODO
- std::vector<Dof> tkeys;
- this->getKeys(ele,tkeys);
- int nkeys = this->getNumKeys(ele);
- int ncomp = _comp.size();
- int nbver = nkeys/ncomp-1; // because GC is in the nkeys
- for(int i=0;i<ncomp; i++)
- keys.push_back(tkeys[ind+i*nbver]);
- this->getKeysOfGravityCenter(keys);
- }
- virtual int getNumKeys(MElement *ele, int ind){
- return 2*_comp.size();
- }
+ virtual void getKeys(MElement *ele,std::vector<Dof> &keys)=0;
+ virtual int getNumKeys(MElement *ele)=0;
+ virtual void getKeysOfGravityCenter(std::vector<Dof> &keys)=0;
+ virtual int getNumKeysOfGravityCenter()=0;
+ virtual void getKeys(MElement *ele, const int ind, std::vector<Dof> &keys)=0;
+ virtual int getNumKeys(MElement *ele, int ind)=0; // number key in one component + number of key of GC
};
#endif // CONTACTSPACE_H_
diff --git a/NonLinearSolver/contact/contactTerms.h b/NonLinearSolver/contact/contactTerms.h
index d0b652d20f..fbf61943dd 100644
--- a/NonLinearSolver/contact/contactTerms.h
+++ b/NonLinearSolver/contact/contactTerms.h
@@ -50,7 +50,7 @@ template <class T2=double>class contactBilinearTermBase : public BilinearTermBa
template<class T2=double> class rigidContactLinearTermBase : public contactLinearTermBase<T2>{
protected:
const unknownField *_ufield;
- rigidContactSpace *_spc;
+ rigidContactSpaceBase *_spc;
// data for get function (Allocated once)
mutable SVector3 _lastNormalContact; // allow to use the normal in the 2 get functions not very elegant
mutable int nbdofs,nbdofsGC, nbvertex,nbcomp,nbcomptimenbvertex;
@@ -60,7 +60,7 @@ template<class T2=double> class rigidContactLinearTermBase : public contactLinea
mutable double mvalue[6];
mutable MVertex *ver;
public:
- rigidContactLinearTermBase(rigidContactSpace *sp,
+ rigidContactLinearTermBase(rigidContactSpaceBase *sp,
const unknownField *ufield) : contactLinearTermBase<T2>(),
_ufield(ufield), _spc(sp), _lastNormalContact(0.,0.,0.),
nbdofs(0), nbdofsGC(0), nbvertex(0), nbcomp(0), nbcomptimenbvertex(0)
diff --git a/NonLinearSolver/contact/rigidCylinderContactTerms.cpp b/NonLinearSolver/contact/rigidCylinderContactTerms.cpp
index ab732190a1..ed6600e418 100644
--- a/NonLinearSolver/contact/rigidCylinderContactTerms.cpp
+++ b/NonLinearSolver/contact/rigidCylinderContactTerms.cpp
@@ -11,8 +11,8 @@
#include "rigidCylinderContactTerms.h"
#include "contactDomain.h"
#include "unknownField.h"
-massRigidCylinder::massRigidCylinder(rigidCylinderContact *cdom,
- rigidContactSpace *spc) : _length(cdom->getLength()),_radius(cdom->getRadius()),
+massRigidCylinder::massRigidCylinder(rigidCylinderContactDomain *cdom,
+ rigidContactSpaceBase *spc) : _length(cdom->getLength()),_radius(cdom->getRadius()),
_thick(cdom->getThickness()),_rho(cdom->getDensity()),
_spc(spc),nbdofs(0.), masse(0.),Rint(0.){}
@@ -32,8 +32,8 @@ void massRigidCylinder::get(MElement *ele,int npts,IntPt *GP,fullMatrix<double>
m(i,i) = masse;
}
- forceRigidCylinderContact::forceRigidCylinderContact(const rigidCylinderContact *cdom,
- rigidContactSpace *sp,const double thick,
+ forceRigidCylinderContact::forceRigidCylinderContact(const rigidCylinderContactDomain *cdom,
+ rigidContactSpaceBase *sp,const double thick,
const unknownField *ufield) : rigidContactLinearTermBase(sp,ufield), _cdom(cdom),
_vergc(cdom->getGC()),
_axisDirector(cdom->getAxisDirector()),
diff --git a/NonLinearSolver/contact/rigidCylinderContactTerms.h b/NonLinearSolver/contact/rigidCylinderContactTerms.h
index de8a87eec7..ba31dd9248 100644
--- a/NonLinearSolver/contact/rigidCylinderContactTerms.h
+++ b/NonLinearSolver/contact/rigidCylinderContactTerms.h
@@ -13,18 +13,18 @@
#define RIGIDCONTACTCYLINDER_H_
#include "contactTerms.h"
#include "contactFunctionSpace.h"
-class rigidCylinderContact;
+class rigidCylinderContactDomain;
class massRigidCylinder : public BilinearTermBase{
protected:
double _length, _radius, _thick, _rho;
- rigidContactSpace *_spc;
+ rigidContactSpaceBase *_spc;
// Data for get function (Allocated once)
mutable int nbdofs;
mutable double masse;
mutable double Rint;
public:
- massRigidCylinder(rigidCylinderContact *cdom,rigidContactSpace *spc);
- massRigidCylinder(rigidContactSpace *spc,double leng, double radius,
+ massRigidCylinder(rigidCylinderContactDomain *cdom,rigidContactSpaceBase *spc);
+ massRigidCylinder(rigidContactSpaceBase *spc,double leng, double radius,
double thick,double rho) : _spc(spc), _length(leng), _radius(radius),_thick(thick), _rho(rho){}
// Arguments are useless but it's allow to use classical assemble function
void get(MElement *ele, int npts, IntPt *GP, fullMatrix<double> &m) const;
@@ -40,7 +40,7 @@ class massRigidCylinder : public BilinearTermBase{
class forceRigidCylinderContact : public rigidContactLinearTermBase<double>{
protected:
- const rigidCylinderContact *_cdom;
+ const rigidCylinderContactDomain *_cdom;
double _fdgfactor, _facGC;
const MVertex *_vergc;
const SVector3* _axisDirector;
@@ -54,8 +54,8 @@ class forceRigidCylinderContact : public rigidContactLinearTermBase<double>{
mutable double d,penetration,penpen;
mutable SVector3 dirAC;
public:
- forceRigidCylinderContact(const rigidCylinderContact *cdom,
- rigidContactSpace *sp,const double thick,
+ forceRigidCylinderContact(const rigidCylinderContactDomain *cdom,
+ rigidContactSpaceBase *sp,const double thick,
const unknownField *ufield);
~forceRigidCylinderContact(){}
// void get(MElement *ele,int npts, IntPt *GP,fullVector<double> &m);
@@ -79,7 +79,7 @@ class stiffnessRigidCylinderContact : public contactBilinearTermBase<double>{
void get(MElement *ele, const int verIndex, fullMatrix<double> &m) const;
protected:
const unknownField *_ufield;
- rigidContactSpace *_spc;
+ rigidContactSpaceBase *_spc;
const double _perturbation; // numerical perturbation
const double _doublepertexpm1;
// Data for get function (Allocated Once)
@@ -91,7 +91,7 @@ class stiffnessRigidCylinderContact : public contactBilinearTermBase<double>{
mutable std::vector<double> disp;
mutable double pertdisp[6];
public:
- stiffnessRigidCylinderContact(rigidContactSpace *space, contactLinearTermBase<double> *lterm, const unknownField *ufield) : contactBilinearTermBase(lterm),
+ stiffnessRigidCylinderContact(rigidContactSpaceBase *space, contactLinearTermBase<double> *lterm, const unknownField *ufield) : contactBilinearTermBase(lterm),
_ufield(ufield),
_perturbation(1.e-10),
_doublepertexpm1(1./(2.e-10)),
diff --git a/NonLinearSolver/internalPoints/ipField.h b/NonLinearSolver/internalPoints/ipField.h
index a5a427f36a..86866f7bfd 100644
--- a/NonLinearSolver/internalPoints/ipField.h
+++ b/NonLinearSolver/internalPoints/ipField.h
@@ -12,7 +12,6 @@
# ifndef _IPFIELD_H_
# define _IPFIELD_H_
#include<vector>
-#include"Dof3IntType.h"
#include"quadratureRules.h"
#include "unknownField.h"
#include "elementField.h"
@@ -121,7 +120,7 @@ class IPField : public elementField {
template<class T1,class T2> void getIPv(const MInterfaceElement *iele, const T1** vipv_m, const T2** vipv_p,
const IPStateBase::whichState ws=IPStateBase::current) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size()/2;
for(int i=0;i<npts;i++){
@@ -134,7 +133,7 @@ class IPField : public elementField {
template<class T1,class T2> void getIPv(const MInterfaceElement *iele, const T1** vipv_m, const T1** vipvprev_m,
const T2** vipv_p, const T2** vipvprev_p) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size()/2;
for(int i=0;i<npts;i++){
@@ -149,7 +148,7 @@ class IPField : public elementField {
template<class T1,class T2> void getIPv(const MInterfaceElement *iele, T1** vipv_m, T2** vipv_p,
const IPStateBase::whichState ws=IPStateBase::current) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size()/2;
for(int i=0;i<npts;i++){
@@ -162,7 +161,7 @@ class IPField : public elementField {
template<class T1,class T2> void getIPv(const MInterfaceElement *iele, T1** vipv_m, T2** vipvprev_m,
T2** vipv_p, T2** vipvprev_p) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size()/2;
for(int i=0;i<npts;i++){
@@ -179,7 +178,7 @@ class IPField : public elementField {
template<class T1> void getIPv(const MInterfaceElement *iele, const T1** vipv_m,
const IPStateBase::whichState ws=IPStateBase::current) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size();
for(int i=0;i<npts;i++){
@@ -189,7 +188,7 @@ class IPField : public elementField {
}
template<class T1> void getIPv(const MInterfaceElement *iele, const T1** vipv_m, const T1** vipvprev_m) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
// SAME NUMBER OF GAUSS POINTS ON BOTH SIDES
int npts = vips->size();
for(int i=0;i<npts;i++){
@@ -202,7 +201,7 @@ class IPField : public elementField {
template<class T1> void getIPv(const MInterfaceElement *iele, T1** vipv_m,
const IPStateBase::whichState ws=IPStateBase::current) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
int npts = vips->size();
for(int i=0;i<npts;i++){
IPStateBase *ipsm = (*vips)[i];
@@ -211,7 +210,7 @@ class IPField : public elementField {
}
template<class T1> void getIPv(const MInterfaceElement *iele, T1** vipv_m, T1** vipvprev_m) const
{
- AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNum());
+ AllIPState::ipstateElementContainer *vips = _AIPS->getIPstate(iele->getNumber());
int npts = vips->size();
for(int i=0;i<npts;i++){
IPStateBase *ipsm = (*vips)[i];
diff --git a/NonLinearSolver/internalPoints/ipstate.cpp b/NonLinearSolver/internalPoints/ipstate.cpp
index 0d7ab3c965..f2197c2253 100644
--- a/NonLinearSolver/internalPoints/ipstate.cpp
+++ b/NonLinearSolver/internalPoints/ipstate.cpp
@@ -77,7 +77,7 @@ AllIPState::AllIPState(GModel *pModel, std::vector<partDomain*> &vdom)
for(int i=npts_inter/2;i<npts_inter;i++){
mlawPlus->createIPState(tp[i],&state,ele,iele->getElem(1)->getNumVertices());
}
- _mapall.insert(ipstatePairType(iele->getNum(),tp));
+ _mapall.insert(ipstatePairType(iele->getNumber(),tp));
}
// Virtual interface element (no duplication)
materialLaw *mlaw = dgdom->getMaterialLaw();
@@ -89,7 +89,7 @@ AllIPState::AllIPState(GModel *pModel, std::vector<partDomain*> &vdom)
for(int i=0;i<npts_inter;i++){
mlaw->createIPState(tp[i],&state,ele,iele->getElem(0)->getNumVertices());
}
- _mapall.insert(ipstatePairType(iele->getNum(),tp));
+ _mapall.insert(ipstatePairType(iele->getNumber(),tp));
}
}
// bulk element
diff --git a/NonLinearSolver/nlsolver/CMakeLists.txt b/NonLinearSolver/nlsolver/CMakeLists.txt
index ef986dbe12..ec91530b0d 100644
--- a/NonLinearSolver/nlsolver/CMakeLists.txt
+++ b/NonLinearSolver/nlsolver/CMakeLists.txt
@@ -9,7 +9,6 @@ set(SRC
unknownField.cpp
unknownStaticField.cpp
# inline
- Dof3IntType.h
explicitDofManager.h
explicitHulbertChung.h
nlsolAlgorithms.h
diff --git a/NonLinearSolver/nlsolver/Dof3IntType.h b/NonLinearSolver/nlsolver/Dof3IntType.h
deleted file mode 100644
index 0c7aba525e..0000000000
--- a/NonLinearSolver/nlsolver/Dof3IntType.h
+++ /dev/null
@@ -1,32 +0,0 @@
-//
-// Description: Derivate class of Dof used for FullDg formulation where the type includes 3 informations (field,local vertex number, dof)
-//
-//
-// Author: <Gauthier BECKER>, (C) 2010
-//
-// Copyright: See COPYING file that comes with this distribution
-//
-//
-#ifndef DGC0DOF_H_
-#define DGC0DOF_H_
-#include "dofManager.h"
-// Two functions are added to dof to define type with three int
-// Use by high level so impossible to let this in the projects ?
-class Dof3IntType : public Dof{
- public :
- Dof3IntType(long int ent, int type) : Dof(ent,type){}
- ~Dof3IntType(){}
- inline static int createTypeWithThreeInts(int comp,int field, int num=0){return 100000*field+100*num+comp;}
- inline static void getThreeIntsFromType(int t, int &comp, int &field, int &num){
- field = t/100000;
- int i1=t%100000;
- num = i1/100;
- comp= i1%100;
- }
- // This function is used for Cg/Dg case (Allow to call the same function for the Two formulations)
- inline static void getThreeIntsFromType(int t, int &comp, int &field){
- field = t/100000;
- comp =t%100000;
- }
-};
-#endif // Dof3IntType
diff --git a/NonLinearSolver/nlsolver/nlsolAlgorithms.h b/NonLinearSolver/nlsolver/nlsolAlgorithms.h
index 6f922cf652..d74b712fe9 100644
--- a/NonLinearSolver/nlsolver/nlsolAlgorithms.h
+++ b/NonLinearSolver/nlsolver/nlsolAlgorithms.h
@@ -82,7 +82,7 @@ template<class Iterator,class Assembler> void AssembleMass(BilinearTermBase *ter
// Assemble mass for a rigid contact space. Only ddl of GC
template<class Assembler>
-void AssembleMass(BilinearTermBase *term, rigidContactSpace *space,Assembler *assembler){
+void AssembleMass(BilinearTermBase *term, rigidContactSpaceBase *space,Assembler *assembler){
fullMatrix<typename Assembler::dataMat> localMatrix;
explicitHCDofManager<typename Assembler::dataVec>* expAss = dynamic_cast<explicitHCDofManager<typename Assembler::dataVec>*>(assembler);
std::vector<Dof> R;
@@ -146,7 +146,7 @@ template<class Iterator,class Assembler> void FixNodalDofs(FunctionSpaceBase *sp
}
template<class Assembler>
-void FixNodalDofs(rigidContactSpace *space,simpleFunction<typename Assembler::dataVec> &fct,FilterDof &filter, Assembler &assembler){
+void FixNodalDofs(rigidContactSpaceBase *space,simpleFunction<typename Assembler::dataVec> &fct,FilterDof &filter, Assembler &assembler){
std::vector<Dof> R;
space->getKeysOfGravityCenter(R);
for(int i=0;i<R.size();i++){
@@ -214,7 +214,7 @@ template<class Iterator,class Assembler> void SetInitialDofs(FunctionSpaceBase *
// Function Numbering Dof for rigid contact (Create Three Dofs for GC of rigid bodies)
template<class Assembler>
-void NumberDofs(rigidContactSpace &space, Assembler &assembler){
+void NumberDofs(rigidContactSpaceBase &space, Assembler &assembler){
// get Dofs of GC
std::vector<Dof> R;
space.getKeysOfGravityCenter(R);
diff --git a/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp b/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
index 3e7fdce439..29a9321758 100644
--- a/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
+++ b/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
@@ -91,7 +91,7 @@ void nonLinearMechSolver::initContactInteraction(){
MElement *ele = *ite;
if(ele == elementcont){
flagfind =true;
- cdom->setDomain(dom);
+ cdom->setDomainAndFunctionSpace(dom);
// spaceManager->addSpace(cdom);
break;
}
@@ -107,7 +107,7 @@ void nonLinearMechSolver::initContactInteraction(){
for(contactContainer::iterator itcdom = _allContact.begin(); itcdom!=_allContact.end(); ++itcdom){
contactDomain *cdom = *itcdom;
if(cdom->getPhys() == diri._tag){
- rigidContactSpace *sp = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
+ rigidContactSpaceBase *sp = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
diri.setSpace(sp);
break;
}
@@ -342,7 +342,7 @@ void nonLinearMechSolver::numberDofs(linearSystem<double> *lsys){
for(contactContainer::iterator it = _allContact.begin(); it!=_allContact.end(); ++it){
contactDomain* cdom = *it;
if(cdom->isRigid()){
- rigidContactSpace *rcspace = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
+ rigidContactSpaceBase *rcspace = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
NumberDofs(*rcspace,*pAssembler);
}
}
@@ -566,13 +566,8 @@ void nonLinearMechSolver::initTerms(unknownField *ufield, IPField *ipf){
// contact domain
for(contactContainer::iterator it = _allContact.begin(); it!=_allContact.end(); ++it){
contactDomain *cdom = *it;
- if(cdom->getContactType() == contactDomain::rigidCylinder){
- rigidCylinderContact *rcc = dynamic_cast<rigidCylinderContact*>(cdom);
- rcc->initializeTerms(ufield);
- }
- else Msg::Error("Terms are not defined for contact type %d",cdom->getContactType());
+ cdom->initializeTerms(ufield);
}
-
}
void nonLinearMechSolver::solveStaticLinar()
@@ -1017,9 +1012,9 @@ void nonLinearMechSolver::archivingNode(const int numphys, const int comp, initi
std::vector<MVertex*> vv;
pModel->getMeshVerticesForPhysicalGroup(0,numphys,vv);
- if(vv.size() !=0){
- std::pair< Dof,initialCondition::whichCondition > pai(Dof(vv[0]->getNum(),Dof3IntType::createTypeWithThreeInts(comp,_tag)),wc);
- anoded.push_back(pai);
+ if(vv.size() == 1){
+// std::pair< Dof,initialCondition::whichCondition > pai(Dof(vv[0]->getNum(),Dof3IntType::createTypeWithThreeInts(comp,_tag)),wc);
+ anoded.push_back(archiveDispNode(vv[0]->getNum(),comp,_tag,wc));
// remove old file (linux only ??)
std::ostringstream oss;
oss << vv[0]->getNum(); // Change this
@@ -1042,7 +1037,10 @@ void nonLinearMechSolver::archivingNode(const int numphys, const int comp, initi
system(rfname.c_str());
}
else{
- Msg::Warning("No physical group %d So it is impossible to archive nodal displacement",numphys);
+ if(vv.size()==0)
+ Msg::Warning("No physical group %d So it is impossible to archive nodal displacement",numphys);
+ else
+ Msg::Warning("More than one node in physical group impssible to archive for now (no loop)",numphys);
}
}
@@ -1504,7 +1502,7 @@ void nonLinearMechSolver::solveExplicit(){
for(nonLinearMechSolver::contactContainer::iterator itc = _allContact.begin(); itc!= _allContact.end(); ++itc){
contactDomain *cdom = *itc;
if(cdom->isRigid()){
- rigidContactSpace *rcsp = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
+ rigidContactSpaceBase *rcsp = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
AssembleMass(cdom->getMassTerm(),rcsp,pAssembler);
}
}
diff --git a/NonLinearSolver/nlsolver/nonLinearMechSolver.h b/NonLinearSolver/nlsolver/nonLinearMechSolver.h
index 50f9b75e19..3ee0deea76 100644
--- a/NonLinearSolver/nlsolver/nonLinearMechSolver.h
+++ b/NonLinearSolver/nlsolver/nonLinearMechSolver.h
@@ -42,6 +42,15 @@ struct archiveForce{
~archiveForce(){}
};
+struct archiveDispNode{
+ int _vernum;
+ int _comp;
+ initialCondition::whichCondition _wc;
+ int _tag;
+ archiveDispNode(const int vernum, const int comp, const int tag,initialCondition::whichCondition wc) : _vernum(vernum), _comp(comp), _wc(wc), _tag(tag){}
+ ~archiveDispNode(){}
+};
+
struct archiveRigidContactDisp{
int _physmaster;
int _comp;
@@ -124,7 +133,7 @@ class nonLinearMechSolver
std::map<int,double> aefvalue;
std::vector<fextPrescribedDisp> _allaef;
// std vector to archive a node displacement
- std::vector<std::pair<Dof,initialCondition::whichCondition> > anoded;
+ std::vector<archiveDispNode> anoded;
// std::vector to archive a force (info musy be store before creation because functionSpace are not initialize)
std::vector<archiveForce> vaf;
// std vector to archive ipvariable
diff --git a/NonLinearSolver/nlsolver/unknownDynamicField.cpp b/NonLinearSolver/nlsolver/unknownDynamicField.cpp
index 38fe57ceb3..967dfdf5a2 100644
--- a/NonLinearSolver/nlsolver/unknownDynamicField.cpp
+++ b/NonLinearSolver/nlsolver/unknownDynamicField.cpp
@@ -13,8 +13,9 @@
#include "nlsolAlgorithms.h"
unknownDynamicField::unknownDynamicField(dofManager<double> *pas, std::vector<partDomain*> &vdom,
nonLinearMechSolver::contactContainer *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
- std::vector<archiveRigidContactDisp> &contactarch, const bool view_, const std::string filen) : unknownField(pas,vdom,acontact,nc,archiving,contactarch,view_,filen)
+ const int nc, std::vector<archiveDispNode> &archiving,
+ std::vector<archiveRigidContactDisp> &contactarch, const bool view_, const std::string filen) : unknownField(pas,vdom,acontact,nc,archiving,
+ contactarch,view_,filen)
{
// to avoid constant dynamic cast
dynassembler = dynamic_cast<explicitHCDofManager<double>*>(pAssembler);
@@ -149,8 +150,7 @@ void unknownDynamicField::archiving(const double time){
oss << it->nodenum;
std::string s = oss.str();
// component of displacement
- int field,comp,num;
- Dof3IntType::getThreeIntsFromType(it->D.getType(),comp,field,num);
+ int comp = it->_comp;
oss.str("");
oss << comp;
std::string s2 = oss.str();
diff --git a/NonLinearSolver/nlsolver/unknownDynamicField.h b/NonLinearSolver/nlsolver/unknownDynamicField.h
index e42eae29a7..7314d996fd 100644
--- a/NonLinearSolver/nlsolver/unknownDynamicField.h
+++ b/NonLinearSolver/nlsolver/unknownDynamicField.h
@@ -21,7 +21,7 @@ class unknownDynamicField : public unknownField{
public:
unknownDynamicField(dofManager<double> *pas, std::vector<partDomain*> &vdom,
nonLinearMechSolver::contactContainer *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
+ const int nc, std::vector<archiveDispNode> &archiving,
std::vector<archiveRigidContactDisp> &contactarch, const bool view_=true, const std::string filen="disp.msh");
virtual void update();
diff --git a/NonLinearSolver/nlsolver/unknownField.cpp b/NonLinearSolver/nlsolver/unknownField.cpp
index a08cb504e2..027aa8bade 100644
--- a/NonLinearSolver/nlsolver/unknownField.cpp
+++ b/NonLinearSolver/nlsolver/unknownField.cpp
@@ -17,7 +17,7 @@
// constructor
unknownField::unknownField(dofManager<double> *pas, std::vector<partDomain*> &vdom, nonLinearMechSolver::contactContainer *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
+ const int nc, std::vector<archiveDispNode> &archiving,
std::vector<archiveRigidContactDisp> &contactarch, const bool view_, const std::string filen): pAssembler(pas),
domainVector(&vdom),
_allContact(acontact),
@@ -43,15 +43,12 @@ unknownField::unknownField(dofManager<double> *pas, std::vector<partDomain*> &vd
sp->getKeys(&ele,R);
// make dof for archiving node
for(int i=0;i<archiving.size();i++){
- if( (alfind[i] == false) and (ver->getNum() == archiving[i].first.getEntity())){
- // get the comp of the archiving node
- int comp,field,num;
- Dof3IntType::getThreeIntsFromType(archiving[i].first.getType(),comp,field,num);
- FilterDof* filter = dom->createFilterDof(comp);
+ if( (alfind[i] == false) and (ver->getNum() == archiving[i]._vernum)){
+ FilterDof* filter = dom->createFilterDof(archiving[i]._comp);
for(int j=0;j<R.size();j++){
if(filter->operator()(R[j])){
alfind[i] = true;
- varch.push_back(archiveNode(R[j],ver->getNum(),archiving[i].second));
+ varch.push_back(archiveNode(R[j],ver->getNum(),archiving[i]._comp,archiving[i]._wc));
}
}
}
@@ -69,12 +66,11 @@ unknownField::unknownField(dofManager<double> *pas, std::vector<partDomain*> &vd
for(int j=0;j<nbvertex;j++){
MVertex *ver = ele->getVertex(j);
for(int i=0; i<archiving.size(); i++){
- if((alfind[i] == false) and(ver->getNum() == archiving[i].first.getEntity())){
+ if((alfind[i] == false) and(ver->getNum() == archiving[i]._vernum)){
// get the comp of the archiving node
- int comp,field,num;
- Dof3IntType::getThreeIntsFromType(archiving[i].first.getType(),comp,field,num);
alfind[i] = true;
- varch.push_back(archiveNode(R[j+comp*nbvertex],ver->getNum(),archiving[i].second));
+ int comp = archiving[i]._comp;
+ varch.push_back(archiveNode(R[j+comp*nbvertex],ver->getNum(),comp,archiving[i]._wc));
break;
}
}
@@ -92,7 +88,7 @@ unknownField::unknownField(dofManager<double> *pas, std::vector<partDomain*> &vd
int pMaster = cdom->getTag();
for(int i=0;i<contactarch.size(); i++){
if(contactarch[i]._physmaster == pMaster){
- rigidContactSpace *sp = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
+ rigidContactSpaceBase *sp = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
std::vector<Dof> R;
sp->getKeysOfGravityCenter(R);
for(int j=0;j<R.size(); j++)
@@ -104,3 +100,81 @@ unknownField::unknownField(dofManager<double> *pas, std::vector<partDomain*> &vd
}
this->setTotElem(totelem);
}
+
+// Add contact archiving
+void unknownField::buildView(std::vector<partDomain*> &vdom,const double time,
+ const int nstep, const std::string &valuename, const int cc, const bool binary){
+ if(view){
+ // test file size (and create a new one if needed)
+ uint32_t fileSize;
+ FILE *fp = fopen(fileName.c_str(), "r");
+ if(!fp){
+ Msg::Error("Unable to open file '%s'", fileName.c_str());
+ return;
+ }
+ fseek (fp, 0, SEEK_END);
+ fileSize = (uint32_t) (ftell(fp));
+ if(fileSize > fmaxsize) this->updateFileName();
+ fclose(fp);
+ fp = fopen(fileName.c_str(), "a");
+
+ // compute the number of element
+ if(binary) Msg::Warning("Binary write not implemented yet");
+ fprintf(fp, "$%s\n",dataname.c_str());
+ fprintf(fp, "1\n\"%s\"\n",valuename.c_str());
+ fprintf(fp, "1\n%.16g\n", time);
+ fprintf(fp, "3\n%d\n%d\n%d\n", nstep, numcomp, totelem);
+ std::vector<double> fieldData;
+ if(type == ElementNodeData){
+ for(std::vector<partDomain*>::iterator itdom=vdom.begin(); itdom!=vdom.end(); ++itdom){
+ partDomain *dom = *itdom;
+ for (groupOfElements::elementContainer::const_iterator it = dom->g->begin(); it != dom->g->end(); ++it){
+ MElement *ele = *it;
+ int numv = ele->getNumVertices();
+ this->get(dom,ele,fieldData,cc);
+ fprintf(fp, "%d %d",ele->getNum(),numv);
+ for(int i=0;i<numv;i++)
+ for(int j=0;j<numcomp;j++)
+ fprintf(fp, " %.16g",fieldData[i+j*numv]);
+ fprintf(fp,"\n");
+ fieldData.clear();
+ }
+ }
+ // loop on contact domain
+ for(nonLinearMechSolver::contactContainer::iterator it=_allContact->begin(); it!=_allContact->end(); ++it){
+ contactDomain *cdom = *it;
+ if(cdom->isRigid()){
+ for (groupOfElements::elementContainer::const_iterator it = cdom->gMaster->begin(); it != cdom->gMaster->end(); ++it){
+ MElement *ele = *it;
+ int numv = ele->getNumVertices();
+ rigidContactSpaceBase *spgc = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
+ std::vector<Dof> R;
+ spgc->getKeysOfGravityCenter(R);
+ this->get(R,fieldData);
+ fprintf(fp, "%d %d",ele->getNum(),numv);
+ for(int i=0;i<numv;i++)
+ for(int j=0;j<numcomp;j++)
+ fprintf(fp, " %.16g",fieldData[j]);
+ fprintf(fp,"\n");
+ fieldData.clear();
+ }
+ }
+ }
+ }
+ else if(type == ElementData){
+ for (unsigned int i = 0; i < vdom.size(); ++i)
+ for (groupOfElements::elementContainer::const_iterator it = vdom[i]->g->begin(); it != vdom[i]->g->end(); ++it){
+ MElement *ele = *it;
+ fieldData.resize(numcomp);
+ this->get(vdom[i],ele,fieldData,cc);
+ fprintf(fp, "%d",ele->getNum());
+ for(int j=0;j<numcomp;j++)
+ fprintf(fp, " %.16g",fieldData[j]);
+ fprintf(fp,"\n");
+ }
+ }
+ fprintf(fp, "$End%s\n",dataname.c_str());
+ fclose(fp);
+ }
+ else Msg::Warning("No displacement view created because the variable view is set to false for this field\n");
+}
diff --git a/NonLinearSolver/nlsolver/unknownField.h b/NonLinearSolver/nlsolver/unknownField.h
index 93c510554d..e15255dc5d 100644
--- a/NonLinearSolver/nlsolver/unknownField.h
+++ b/NonLinearSolver/nlsolver/unknownField.h
@@ -23,11 +23,13 @@
class contactDomain;
class elementField;
struct archiveRigidContactDisp;
+struct archiveDispNode;
struct archiveNode{
long int nodenum;
Dof D;
+ int _comp;
initialCondition::whichCondition wc;
- archiveNode(Dof R, int n, initialCondition::whichCondition wv=initialCondition::position) : nodenum(n), D(R), wc(wv){}
+ archiveNode(Dof R, int n, int comp,initialCondition::whichCondition wv=initialCondition::position) : nodenum(n), D(R), _comp(comp), wc(wv){}
~archiveNode(){};
};
@@ -42,7 +44,7 @@ class unknownField : public elementField{
public:
// update all displacement value
unknownField(dofManager<double> *pas, std::vector<partDomain*> &elas, std::set<contactDomain*> *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
+ const int nc, std::vector<archiveDispNode> &archiving,
std::vector<archiveRigidContactDisp> &contactarch, const bool =true, const std::string="disp.msh");
~unknownField(){}
virtual void update()=0;
@@ -58,7 +60,7 @@ class unknownField : public elementField{
virtual void get(FunctionSpaceBase *sp1,FunctionSpaceBase *sp2, MInterfaceElement *iele,std::vector<double> &udofs)=0;
virtual void archiving(const double time)=0;
virtual void setInitial(const std::vector<Dof> &R,const std::vector<double> &disp)=0;
-// virtual void buildView(std::vector<partDomain*> &vdom,const double time,
-// const int nstep, const std::string &valuename, const int cc,const bool binary);
+ virtual void buildView(std::vector<partDomain*> &vdom,const double time,
+ const int nstep, const std::string &valuename, const int cc=-1,const bool binary=false);
};
#endif // _UNKNOWNFIELD_H_
diff --git a/NonLinearSolver/nlsolver/unknownStaticField.cpp b/NonLinearSolver/nlsolver/unknownStaticField.cpp
index 934d5c95eb..2f3a58cd50 100644
--- a/NonLinearSolver/nlsolver/unknownStaticField.cpp
+++ b/NonLinearSolver/nlsolver/unknownStaticField.cpp
@@ -14,7 +14,7 @@
#include "nonLinearMechSolver.h"
#include "MPoint.h"
unknownStaticField::unknownStaticField(dofManager<double> *pas, std::vector<partDomain*> &vdom, nonLinearMechSolver::contactContainer *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
+ const int nc, std::vector<archiveDispNode> &archiving,
std::vector<archiveRigidContactDisp> &contactarch, const bool view_, const std::string filen
) : unknownField(pas,vdom,acontact,nc,archiving,contactarch,view_,filen)
{
@@ -37,7 +37,7 @@ unknownStaticField::unknownStaticField(dofManager<double> *pas, std::vector<part
contactDomain *cdom = *it;
if(cdom->isRigid()){
std::vector<Dof> R;
- rigidContactSpace *sp = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
+ rigidContactSpaceBase *sp = dynamic_cast<rigidContactSpaceBase*>(cdom->getSpace());
sp->getKeysOfGravityCenter(R);
for(int i=0;i<R.size(); i++)
umap.insert(std::pair<Dof,double>(R[i],0.));
@@ -109,17 +109,6 @@ void unknownStaticField::get(std::vector<Dof> &R, fullVector<double> &disp) cons
}
}
-/*void unknownStaticField:get(MVertex *ver,std::vector<double> &udofs){
- long int ent;
- if(!fullDg){
- ent = ver->getNum();
- udofs = umap.find(ent)->second;
- }
- else{
- Msg::Error("Impossible to get displacement by vertex for full Dg formulation. Work only with MElement");
- }
-}*/
-
void unknownStaticField::get(partDomain *dom, MElement *ele,std::vector<double> &udofs, const int cc){
std::vector<Dof> R;
FunctionSpaceBase *sp = dom->getFunctionSpace();
@@ -160,27 +149,6 @@ void unknownStaticField::get(FunctionSpaceBase *sp1, FunctionSpaceBase *sp2,MInt
}
}
-/*void unknownStaticField:getForPerturbation(FunctionSpaceBase &sp, MInterfaceElement* iele, const bool minus, Dof &D, double pert, std::vector<double> &udofs){
- double eps = -pert;
- std::vector<Dof> R;
- sp.getKeys(iele->getElem(0),R);
- if(!minus){
- this->set(D,eps);
- this->get(R,udofs);
- this->set(D,pert);
- }
- else{
- this->get(R,udofs);
- this->set(D,eps);
- }
- if(!(iele->getElem(0)==iele->getElem(1)))// Virtual interface element
- R.clear();
- sp.getKeys(iele->getElem(1),R);
- this->get(R,udofs);
- if(minus)
- this->set(D,pert);
-}*/
-
void unknownStaticField::updateFixedDof(){
double u;
for(std::vector<Dof>::iterator itD=fixedDof.begin(); itD<fixedDof.end();++itD){
@@ -200,8 +168,7 @@ void unknownStaticField::archiving(const double time){
oss << it->nodenum;
std::string s = oss.str();
// component of displacement
- int field,comp,num;
- Dof3IntType::getThreeIntsFromType(it->D.getType(),comp,field,num);
+ int comp = it->_comp;
oss.str("");
oss << comp;
std::string s2 = oss.str();
@@ -217,82 +184,3 @@ void unknownStaticField::setInitial(const std::vector<Dof> &R, const std::vector
umap[R[i]] = disp[i];
}
}
-
-/*
-void unknownField::buildView(std::vector<partDomain*> &vdom,const double time,
- const int nstep, const std::string &valuename, const int cc=-1, const bool binary=false){
- if(view){
- // test file size (and create a new one if needed)
- uint32_t fileSize;
- FILE *fp = fopen(fileName.c_str(), "r");
- if(!fp){
- Msg::Error("Unable to open file '%s'", fileName.c_str());
- return;
- }
- fseek (fp, 0, SEEK_END);
- fileSize = (uint32_t) (ftell(fp));
- if(fileSize > fmaxsize) this->updateFileName();
- fclose(fp);
- fp = fopen(fileName.c_str(), "a");
-
- // compute the number of element
- if(binary) Msg::Warning("Binary write not implemented yet");
- fprintf(fp, "$%s\n",dataname.c_str());
- fprintf(fp, "1\n\"%s\"\n",valuename.c_str());
- fprintf(fp, "1\n%.16g\n", time);
- fprintf(fp, "3\n%d\n%d\n%d\n", nstep, numcomp, totelem);
- std::vector<double> fieldData;
- if(type == ElementNodeData){
- for(std::vector<partDomain*>::iterator itdom=vdom.begin(); itdom!=vdom.end(); ++itdom){
- partDomain *dom = *itdom;
- for (groupOfElements::elementContainer::const_iterator it = dom->g->begin(); it != dom->g->end(); ++it){
- MElement *ele = *it;
- int numv = ele->getNumVertices();
- this->get(dom,ele,fieldData,cc);
- fprintf(fp, "%d %d",ele->getNum(),numv);
- for(int i=0;i<numv;i++)
- for(int j=0;j<numcomp;j++)
- fprintf(fp, " %.16g",fieldData[i+j*numv]);
- fprintf(fp,"\n");
- fieldData.clear();
- }
- }
- // loop on contact domain
- for(nonLinearMechSolver::contactContainer::iterator it=_allContact->begin(); it!=_allContact->end(); ++it){
- contactDomain *cdom = *it;
- if(cdom->isRigid()){
- for (groupOfElements::elementContainer::const_iterator it = cdom->gMaster->begin(); it != cdom->gMaster->end(); ++it){
- MElement *ele = *it;
- int numv = ele->getNumVertices();
- rigidContactSpace *spgc = dynamic_cast<rigidContactSpace*>(cdom->getSpace());
- std::vector<Dof> R;
- spgc->getKeysOfGravityCenter(R);
- this->get(R,fieldData);
- fprintf(fp, "%d %d",ele->getNum(),numv);
- for(int i=0;i<numv;i++)
- for(int j=0;j<numcomp;j++)
- fprintf(fp, " %.16g",fieldData[j]);
- fprintf(fp,"\n");
- fieldData.clear();
- }
- }
- }
- }
- else if(type == ElementData){
- for (unsigned int i = 0; i < vdom.size(); ++i)
- for (groupOfElements::elementContainer::const_iterator it = vdom[i]->g->begin(); it != vdom[i]->g->end(); ++it){
- MElement *ele = *it;
- fieldData.resize(numcomp);
- this->get(vdom[i],ele,fieldData,cc);
- fprintf(fp, "%d",ele->getNum());
- for(int j=0;j<numcomp;j++)
- fprintf(fp, " %.16g",fieldData[j]);
- fprintf(fp,"\n");
- }
- }
- fprintf(fp, "$End%s\n",dataname.c_str());
- fclose(fp);
- }
- else Msg::Warning("No displacement view created because the variable view is set to false for this field\n");
-}
-*/
diff --git a/NonLinearSolver/nlsolver/unknownStaticField.h b/NonLinearSolver/nlsolver/unknownStaticField.h
index 8ac8d04544..16efbea222 100644
--- a/NonLinearSolver/nlsolver/unknownStaticField.h
+++ b/NonLinearSolver/nlsolver/unknownStaticField.h
@@ -25,7 +25,7 @@ class unknownStaticField : public unknownField{
std::map<Dof,double> umap; // One Entry by Dof allow to manage more than 1 domain
public :
unknownStaticField(dofManager<double> *pas, std::vector<partDomain*> &elas, std::set<contactDomain*> *acontact,
- const int nc, std::vector<std::pair<Dof,initialCondition::whichCondition> > &archiving,
+ const int nc, std::vector<archiveDispNode> &archiving,
std::vector<archiveRigidContactDisp> &contactarch, const bool =true, const std::string="disp.msh"
) ;
// update all displacement value
diff --git a/Post/PView.h b/Post/PView.h
index 28006ed823..6db1291f71 100644
--- a/Post/PView.h
+++ b/Post/PView.h
@@ -53,10 +53,10 @@ class PView{
std::vector<double> &x, std::vector<double> &y);
// construct a new model-based view from a bunch of data
PView(std::string name, std::string type, GModel *model,
- std::map<int, std::vector<double> > &data, double time=0.,
+ std::map<int, std::vector<double> > &data, double time=0.,
int numComp = -1);
// add a new time step to a given model-based view
- void addStep(GModel *model, std::map<int, std::vector<double> > &data,
+ void addStep(GModel *model, std::map<int, std::vector<double> > &data,
double time=0.,int numComp = -1);
// default destructor
@@ -66,8 +66,8 @@ class PView{
void deleteVertexArrays();
// get/set the display options
- PViewOptions *getOptions(){ return _options; }
- void setOptions(PViewOptions *val=0);
+ PViewOptions *getOptions(){ return _options; }
+ void setOptions(PViewOptions *val=0);
// get/set the view data
PViewData *getData(bool useAdaptiveIfAvailable=false);
@@ -100,7 +100,7 @@ class PView{
// find view by name or by number (if timeStep >= 0, return view
// only if it does *not* contain that timestep; if partition >= 0,
// return view only if it does *not* contain that partition)
- static PView *getViewByName(std::string name, int timeStep=-1,
+ static PView *getViewByName(std::string name, int timeStep=-1,
int partition=-1);
static PView *getViewByNum(int num, int timeStep=-1, int partition=-1);
@@ -120,7 +120,7 @@ class PView{
void fillVertexArrays();
// fill a vertex array using a raw stream of bytes
- static void fillVertexArray(ConnectionManager *remote, int length,
+ static void fillVertexArray(ConnectionManager *remote, int length,
const char *data, int swap);
// smoothed normals
@@ -131,11 +131,9 @@ class PView{
// this is the maximum number of nodes of elements we actually *draw*
// (high order elements are always subdivided before drawing)
#define PVIEW_NMAX 8
-
void changeCoordinates(PView *p, int ient, int iele,
- int numNodes, int type, int numComp,
+ int numNodes, int type, int numComp,
double **xyz, double **val);
-bool isElementVisible(PViewOptions *opt, int dim, int numNodes,
+bool isElementVisible(PViewOptions *opt, int dim, int numNodes,
double **xyz);
-
#endif
diff --git a/Post/PViewVertexArrays.cpp b/Post/PViewVertexArrays.cpp
index 90f434f987..f6e06bfc1f 100644
--- a/Post/PViewVertexArrays.cpp
+++ b/Post/PViewVertexArrays.cpp
@@ -21,8 +21,8 @@
#include "Options.h"
#include "StringUtils.h"
-static void saturate(int nb, double **val, double vmin, double vmax,
- int i0=0, int i1=1, int i2=2, int i3=3,
+static void saturate(int nb, double **val, double vmin, double vmax,
+ int i0=0, int i1=1, int i2=2, int i3=3,
int i4=4, int i5=5, int i6=6, int i7=7)
{
int id[8] = {i0, i1, i2, i3, i4, i5, i6, i7};
@@ -34,10 +34,10 @@ static void saturate(int nb, double **val, double vmin, double vmax,
static SVector3 normal3(double **xyz, int i0=0, int i1=1, int i2=2)
{
- SVector3 t1(xyz[i1][0] - xyz[i0][0],
+ SVector3 t1(xyz[i1][0] - xyz[i0][0],
xyz[i1][1] - xyz[i0][1],
xyz[i1][2] - xyz[i0][2]);
- SVector3 t2(xyz[i2][0] - xyz[i0][0],
+ SVector3 t2(xyz[i2][0] - xyz[i0][0],
xyz[i2][1] - xyz[i0][1],
xyz[i2][2] - xyz[i0][2]);
SVector3 n = crossprod(t1, t2);
@@ -103,7 +103,7 @@ static void getLineNormal(PView *p, double x[2], double y[2], double z[2],
}
}
-static bool getExternalValues(PView *p, int index, int ient, int iele,
+static bool getExternalValues(PView *p, int index, int ient, int iele,
int numNodes, int numComp, double **val,
int &numComp2, double **val2)
{
@@ -145,7 +145,7 @@ static bool getExternalValues(PView *p, int index, int ient, int iele,
return false;
}
-static void applyGeneralRaise(PView *p, int numNodes, int numComp,
+static void applyGeneralRaise(PView *p, int numNodes, int numComp,
double **vals, double **xyz)
{
PViewOptions *opt = p->getOptions();
@@ -179,10 +179,10 @@ void changeCoordinates(PView *p, int ient, int iele,
xyz[i][j] = barycenter[j] + opt->explode * (xyz[i][j] - barycenter[j]);
}
- if(opt->transform[0][0] != 1. || opt->transform[0][1] != 0. ||
- opt->transform[0][2] != 0. || opt->transform[1][0] != 0. ||
+ if(opt->transform[0][0] != 1. || opt->transform[0][1] != 0. ||
+ opt->transform[0][2] != 0. || opt->transform[1][0] != 0. ||
opt->transform[1][1] != 1. || opt->transform[1][2] != 0. ||
- opt->transform[2][0] != 0. || opt->transform[2][1] != 0. ||
+ opt->transform[2][0] != 0. || opt->transform[2][1] != 0. ||
opt->transform[2][2] != 1.){
for(int i = 0; i < numNodes; i++) {
double old[3] = {xyz[i][0], xyz[i][1], xyz[i][2]};
@@ -241,7 +241,7 @@ void changeCoordinates(PView *p, int ient, int iele,
double **val2 = new double*[numNodes];
for(int i = 0; i < numNodes; i++)
val2[i] = new double[9];
- getExternalValues(p, opt->viewIndexForGenRaise, ient, iele, numNodes,
+ getExternalValues(p, opt->viewIndexForGenRaise, ient, iele, numNodes,
numComp, val, numComp2, val2);
applyGeneralRaise(p, numNodes, numComp2, val2, xyz);
for(int i = 0; i < numNodes; i++)
@@ -252,9 +252,9 @@ void changeCoordinates(PView *p, int ient, int iele,
static double evalClipPlane(int clip, double x, double y, double z)
{
- return CTX::instance()->clipPlane[clip][0] * x +
- CTX::instance()->clipPlane[clip][1] * y +
- CTX::instance()->clipPlane[clip][2] * z +
+ return CTX::instance()->clipPlane[clip][0] * x +
+ CTX::instance()->clipPlane[clip][1] * y +
+ CTX::instance()->clipPlane[clip][2] * z +
CTX::instance()->clipPlane[clip][3];
}
@@ -268,7 +268,7 @@ static double intersectClipPlane(int clip, int numNodes, double **xyz)
return val;
}
-bool isElementVisible(PViewOptions *opt, int dim, int numNodes,
+bool isElementVisible(PViewOptions *opt, int dim, int numNodes,
double **xyz)
{
if(!CTX::instance()->clipWholeElements) return true;
@@ -293,7 +293,7 @@ bool isElementVisible(PViewOptions *opt, int dim, int numNodes,
return !hidden;
}
-static void addOutlinePoint(PView *p, double **xyz, unsigned int color,
+static void addOutlinePoint(PView *p, double **xyz, unsigned int color,
bool pre, int i0=0)
{
if(pre) return;
@@ -301,8 +301,8 @@ static void addOutlinePoint(PView *p, double **xyz, unsigned int color,
p->va_points->add(&xyz[i0][0], &xyz[i0][1], &xyz[i0][2], &n, &color, 0, true);
}
-static void addScalarPoint(PView *p, double **xyz,
- double **val, bool pre, int i0=0,
+static void addScalarPoint(PView *p, double **xyz,
+ double **val, bool pre, int i0=0,
bool unique=false)
{
if(pre) return;
@@ -313,8 +313,8 @@ static void addScalarPoint(PView *p, double **xyz,
if(opt->saturateValues) saturate(1, val, vmin, vmax, i0);
if(val[i0][0] >= vmin && val[i0][0] <= vmax){
- unsigned int col = opt->getColor(val[i0][0], vmin, vmax, false,
- (opt->intervalsType == PViewOptions::Discrete) ?
+ unsigned int col = opt->getColor(val[i0][0], vmin, vmax, false,
+ (opt->intervalsType == PViewOptions::Discrete) ?
opt->nbIso : -1);
SVector3 n = getPointNormal(p, val[i0][0]);
p->va_points->add(&xyz[i0][0], &xyz[i0][1], &xyz[i0][2], &n, &col, 0, unique);
@@ -330,7 +330,7 @@ static void addOutlineLine(PView *p, double **xyz, unsigned int color,
unsigned int col[2];
double x[2], y[2], z[2];
for(int i = 0; i < 2; i++){
- x[i] = xyz[in[i]][0]; y[i] = xyz[in[i]][1]; z[i] = xyz[in[i]][2];
+ x[i] = xyz[in[i]][0]; y[i] = xyz[in[i]][1]; z[i] = xyz[in[i]][2];
col[i] = color;
}
SVector3 n[2];
@@ -338,7 +338,7 @@ static void addOutlineLine(PView *p, double **xyz, unsigned int color,
p->va_lines->add(x, y, z, n, col, 0, true);
}
-static void addScalarLine(PView *p, double **xyz,
+static void addScalarLine(PView *p, double **xyz,
double **val, bool pre, int i0=0, int i1=1,
bool unique=false)
{
@@ -401,7 +401,7 @@ static void addScalarLine(PView *p, double **xyz,
if(vmin == vmax) break;
}
}
-
+
if(opt->intervalsType == PViewOptions::Iso){
for(int k = 0; k < opt->nbIso; k++) {
if(vmin == vmax) k = opt->nbIso / 2;
@@ -537,7 +537,7 @@ static void addScalarTriangle(PView *p, double **xyz, double **val,
if(vmin == vmax) break;
}
}
-
+
if(opt->intervalsType == PViewOptions::Iso){
for(int k = 0; k < opt->nbIso; k++) {
if(vmin == vmax) k = opt->nbIso / 2;
@@ -563,8 +563,8 @@ static void addScalarTriangle(PView *p, double **xyz, double **val,
}
}
-static void addOutlineQuadrangle(PView *p, double **xyz,
- unsigned int color, bool pre, int i0=0, int i1=1,
+static void addOutlineQuadrangle(PView *p, double **xyz,
+ unsigned int color, bool pre, int i0=0, int i1=1,
int i2=2, int i3=3)
{
PViewOptions *opt = p->getOptions();
@@ -590,8 +590,8 @@ static void addOutlineQuadrangle(PView *p, double **xyz,
}
}
-static void addScalarQuadrangle(PView *p, double **xyz,
- double **val, bool pre, int i0=0,
+static void addScalarQuadrangle(PView *p, double **xyz,
+ double **val, bool pre, int i0=0,
int i1=1, int i2=2, int i3=3, bool unique=false)
{
PViewOptions *opt = p->getOptions();
@@ -606,19 +606,19 @@ static void addScalarQuadrangle(PView *p, double **xyz,
opt->boundary++;
return;
}
-
+
for(int i = 0; i < 2; i++)
addScalarTriangle(p, xyz, val, pre, it[i][0], it[i][1], it[i][2], unique);
}
-static void addOutlinePolygon(PView *p, double **xyz,
+static void addOutlinePolygon(PView *p, double **xyz,
unsigned int color, bool pre, int numNodes)
{
for(int i = 0; i < numNodes / 3; i++)
addOutlineTriangle(p, xyz, color, pre, 3*i, 3*i + 1, 3*i + 2);
}
-static void addScalarPolygon(PView *p, double **xyz,
+static void addScalarPolygon(PView *p, double **xyz,
double **val, bool pre, int numNodes)
{
PViewOptions *opt = p->getOptions();
@@ -653,12 +653,12 @@ static void addScalarPolygon(PView *p, double **xyz,
delete verts[i];
return;
}
-
+
for(int i = 0; i < numNodes / 3; i++)
addScalarTriangle(p, xyz, val, pre, 3*i, 3*i+1, 3*i+2);
}
-static void addOutlineTetrahedron(PView *p, double **xyz,
+static void addOutlineTetrahedron(PView *p, double **xyz,
unsigned int color, bool pre)
{
const int it[4][3] = {{0, 2, 1}, {0, 1, 3}, {0, 3, 2}, {3, 1, 2}};
@@ -666,8 +666,8 @@ static void addOutlineTetrahedron(PView *p, double **xyz,
addOutlineTriangle(p, xyz, color, pre, it[i][0], it[i][1], it[i][2]);
}
-static void addScalarTetrahedron(PView *p, double **xyz,
- double **val, bool pre, int i0=0,
+static void addScalarTetrahedron(PView *p, double **xyz,
+ double **val, bool pre, int i0=0,
int i1=1, int i2=2, int i3=3)
{
PViewOptions *opt = p->getOptions();
@@ -722,18 +722,18 @@ static void addScalarTetrahedron(PView *p, double **xyz,
}
}
-static void addOutlineHexahedron(PView *p, double **xyz,
+static void addOutlineHexahedron(PView *p, double **xyz,
unsigned int color, bool pre)
{
const int iq[6][4] = {{0, 3, 2, 1}, {0, 1, 5, 4}, {0, 4, 7, 3},
{1, 2, 6, 5}, {2, 3, 7, 6}, {4, 5, 6, 7}};
for(int i = 0; i < 6; i++)
- addOutlineQuadrangle(p, xyz, color, pre, iq[i][0], iq[i][1],
+ addOutlineQuadrangle(p, xyz, color, pre, iq[i][0], iq[i][1],
iq[i][2], iq[i][3]);
}
-static void addScalarHexahedron(PView *p, double **xyz,
+static void addScalarHexahedron(PView *p, double **xyz,
double **val, bool pre)
{
PViewOptions *opt = p->getOptions();
@@ -750,12 +750,12 @@ static void addScalarHexahedron(PView *p, double **xyz,
opt->boundary++;
return;
}
-
+
for(int i = 0; i < 6; i++)
addScalarTetrahedron(p, xyz, val, pre, is[i][0], is[i][1], is[i][2], is[i][3]);
}
-static void addOutlinePrism(PView *p, double **xyz, unsigned int color,
+static void addOutlinePrism(PView *p, double **xyz, unsigned int color,
bool pre)
{
const int iq[3][4] = {{0, 1, 4, 3}, {0, 3, 5, 2}, {1, 2, 5, 4}};
@@ -766,7 +766,7 @@ static void addOutlinePrism(PView *p, double **xyz, unsigned int color,
for(int i = 0; i < 2; i++)
addOutlineTriangle(p, xyz, color, pre, it[i][0], it[i][1], it[i][2]);
}
-static void addScalarPrism(PView *p, double **xyz,
+static void addScalarPrism(PView *p, double **xyz,
double **val, bool pre)
{
PViewOptions *opt = p->getOptions();
@@ -783,7 +783,7 @@ static void addScalarPrism(PView *p, double **xyz,
opt->boundary++;
return;
}
-
+
for(int i = 0; i < 3; i++)
addScalarTetrahedron(p, xyz, val, pre, is[i][0], is[i][1], is[i][2], is[i][3]);
}
@@ -798,7 +798,7 @@ static void addOutlinePyramid(PView *p, double **xyz,
addOutlineTriangle(p, xyz, color, pre, it[i][0], it[i][1], it[i][2]);
}
-static void addScalarPyramid(PView *p, double **xyz,
+static void addScalarPyramid(PView *p, double **xyz,
double **val, bool pre)
{
PViewOptions *opt = p->getOptions();
@@ -819,7 +819,7 @@ static void addScalarPyramid(PView *p, double **xyz,
addScalarTetrahedron(p, xyz, val, pre, is[i][0], is[i][1], is[i][2], is[i][3]);
}
-static void addOutlinePolyhedron(PView *p, double **xyz,
+static void addOutlinePolyhedron(PView *p, double **xyz,
unsigned int color, bool pre, int numNodes)
{
const int it[4][3] = {{0, 2, 1}, {0, 1, 3}, {0, 3, 2}, {3, 1, 2}};
@@ -847,7 +847,7 @@ static void addOutlinePolyhedron(PView *p, double **xyz,
delete verts[i];
}
-static void addScalarPolyhedron(PView *p, double **xyz,
+static void addScalarPolyhedron(PView *p, double **xyz,
double **val, bool pre, int numNodes)
{
PViewOptions *opt = p->getOptions();
@@ -855,12 +855,12 @@ static void addScalarPolyhedron(PView *p, double **xyz,
if(opt->boundary > 0){
return;
}
-
+
for(int i = 0; i < numNodes / 4; i++)
addScalarTetrahedron(p, xyz, val, pre, 4*i, 4*i + 1, 4*i + 2, 4*i + 3);
}
-static void addOutlineElement(PView *p, int type, double **xyz,
+static void addOutlineElement(PView *p, int type, double **xyz,
bool pre, int numNodes)
{
PViewOptions *opt = p->getOptions();
@@ -895,7 +895,7 @@ static void addScalarElement(PView *p, int type, double **xyz,
}
}
-static void addVectorElement(PView *p, int ient, int iele, int numNodes,
+static void addVectorElement(PView *p, int ient, int iele, int numNodes,
int type, double **xyz,
double **val, bool pre)
{
@@ -907,7 +907,7 @@ static void addVectorElement(PView *p, int ient, int iele, int numNodes,
double **val2 = new double*[numNodes];
for(int i = 0; i < numNodes; i++)
val2[i] = new double[9];
- getExternalValues(p, opt->externalViewIndex, ient, iele, numNodes,
+ getExternalValues(p, opt->externalViewIndex, ient, iele, numNodes,
3, val, numComp2, val2);
if(opt->vectorType == PViewOptions::Displacement){
@@ -938,12 +938,12 @@ static void addVectorElement(PView *p, int ient, int iele, int numNodes,
unsigned int col[2];
double norm[2];
for(int i = 0; i < 2; i++){
- norm[i] = sqrt(dxyz[0][i] * dxyz[0][i] +
- dxyz[1][i] * dxyz[1][i] +
+ norm[i] = sqrt(dxyz[0][i] * dxyz[0][i] +
+ dxyz[1][i] * dxyz[1][i] +
dxyz[2][i] * dxyz[2][i]);
col[i] = opt->getColor(norm[i], opt->tmpMin, opt->tmpMax);
}
- for(int j = 0; j < 3; j++){
+ for(int j = 0; j < 3; j++){
dxyz[j][0] = xyz0[j] + dxyz[j][0] * opt->displacementFactor;
dxyz[j][1] = xyz0[j] + dxyz[j][1] * opt->displacementFactor;
}
@@ -970,7 +970,7 @@ static void addVectorElement(PView *p, int ient, int iele, int numNodes,
double v2 = ComputeScalarRep(numComp2, val2[i]);
if(v2 >= opt->externalMin && v2 <= opt->externalMax){
unsigned int color = opt->getColor(v2, opt->externalMin, opt->externalMax, false,
- (opt->intervalsType == PViewOptions::Discrete) ?
+ (opt->intervalsType == PViewOptions::Discrete) ?
opt->nbIso : -1);
unsigned int col[2] = {color, color};
double dxyz[3][2];
@@ -1003,7 +1003,7 @@ static void addVectorElement(PView *p, int ient, int iele, int numNodes,
if(v2 >= opt->externalMin * (1. - 1.e-15) &&
v2 <= opt->externalMax * (1. + 1.e-15)){
unsigned int color = opt->getColor(v2, opt->externalMin, opt->externalMax, false,
- (opt->intervalsType == PViewOptions::Discrete) ?
+ (opt->intervalsType == PViewOptions::Discrete) ?
opt->nbIso : -1);
unsigned int col[2] = {color, color};
double dxyz[3][2];
@@ -1020,7 +1020,7 @@ static void addVectorElement(PView *p, int ient, int iele, int numNodes,
}
static void addTensorElement(PView *p, int iEnt, int iEle, int numNodes, int type,
- double **xyz, double **val,
+ double **xyz, double **val,
bool pre)
{
PViewOptions *opt = p->getOptions();
@@ -1033,7 +1033,7 @@ static void addTensorElement(PView *p, int iEnt, int iEle, int numNodes, int typ
for(int i = 0; i < numNodes; i++)
val[i][0] = ComputeVonMises(val[i]);
addScalarElement(p, type, xyz, val, pre, numNodes);
- }
+ }
else if (opt->tensorType == PViewOptions::Ellipse ||
opt->tensorType == PViewOptions::Ellipsoid) {
if(opt->glyphLocation == PViewOptions::Vertex){
@@ -1053,7 +1053,7 @@ static void addTensorElement(PView *p, int iEnt, int iEle, int numNodes, int typ
}
double lmax = std::max(S(0), std::max(S(1), S(2)));
unsigned int color = opt->getColor
- (lmax, opt->tmpMin, opt->tmpMax, false,
+ (lmax, opt->tmpMin, opt->tmpMax, false,
(opt->intervalsType == PViewOptions::Discrete) ? opt->nbIso : -1);
unsigned int col[4] = {color, color, color, color};
p->va_ellipses->add(vval[0], vval[1], vval[2], 0, col, 0, false);
@@ -1079,7 +1079,7 @@ static void addTensorElement(PView *p, int iEnt, int iEle, int numNodes, int typ
}
double lmax = std::max(S(0), std::max(S(1), S(2)));
unsigned int color = opt->getColor
- (lmax, opt->tmpMin, opt->tmpMax, false,
+ (lmax, opt->tmpMin, opt->tmpMax, false,
(opt->intervalsType == PViewOptions::Discrete) ? opt->nbIso : -1);
unsigned int col[4] = {color, color, color, color};
p->va_ellipses->add(vval[0], vval[1], vval[2], 0, col, 0, false);
@@ -1131,7 +1131,7 @@ static void addElementsInArrays(PView *p, bool preprocessNormalsOnly)
// use adaptive data if available
PViewData *data = p->getData(true);
PViewOptions *opt = p->getOptions();
-
+
opt->tmpBBox.reset();
int NMAX = PVIEW_NMAX;
@@ -1209,9 +1209,9 @@ static void addElementsInArrays(PView *p, bool preprocessNormalsOnly)
for(int j = 0; j < numNodes; j++)
opt->tmpBBox += SPoint3(xyz[j][0], xyz[j][1], xyz[j][2]);
- if(opt->showElement && !data->useGaussPoints())
+ if(opt->showElement && !data->useGaussPoints())
addOutlineElement(p, type, xyz, preprocessNormalsOnly, numNodes);
-
+
if(opt->intervalsType != PViewOptions::Numeric){
if(data->useGaussPoints()){
for(int j = 0; j < numNodes; j++){
@@ -1256,7 +1256,7 @@ class initPView {
// on Windows/Cygwin
int _estimateIfClipped(PView *p, int num)
{
- if(CTX::instance()->clipWholeElements &&
+ if(CTX::instance()->clipWholeElements &&
CTX::instance()->clipOnlyDrawIntersectingVolume){
PViewOptions *opt = p->getOptions();
for(int clip = 0; clip < 6; clip++){
@@ -1296,10 +1296,10 @@ class initPView {
if(opt->intervalsType == PViewOptions::Iso)
heuristic = (tets + prisms + pyrs + hexas + polyhs) / 10;
else if(opt->intervalsType == PViewOptions::Continuous)
- heuristic = (tris + 2 * quads + 3 * polygs + 6 * tets +
+ heuristic = (tris + 2 * quads + 3 * polygs + 6 * tets +
8 * prisms + 6 * pyrs + 12 * hexas + 10 * polyhs);
else if(opt->intervalsType == PViewOptions::Discrete)
- heuristic = (tris + 2 * quads + 3 * polygs + 6 * tets +
+ heuristic = (tris + 2 * quads + 3 * polygs + 6 * tets +
8 * prisms + 6 * pyrs + 12 * hexas + 10 * polyhs) * 2;
heuristic = _estimateIfClipped(p, heuristic);
return heuristic + 10000;
@@ -1376,11 +1376,11 @@ class initPView {
p->va_vectors->finalize();
p->va_ellipses->finalize();
- Msg::Info("%d vertices in vertex arrays (%g Mb)", p->va_points->getNumVertices() +
- p->va_lines->getNumVertices() + p->va_triangles->getNumVertices() +
+ Msg::Info("%d vertices in vertex arrays (%g Mb)", p->va_points->getNumVertices() +
+ p->va_lines->getNumVertices() + p->va_triangles->getNumVertices() +
p->va_vectors->getNumVertices() + p->va_ellipses->getNumVertices(),
p->va_points->getMemoryInMb() +
- p->va_lines->getMemoryInMb() + p->va_triangles->getMemoryInMb() +
+ p->va_lines->getMemoryInMb() + p->va_triangles->getMemoryInMb() +
p->va_vectors->getMemoryInMb() + p->va_ellipses->getMemoryInMb());
p->setChanged(false);
@@ -1393,7 +1393,7 @@ void PView::fillVertexArrays()
init(this);
}
-void PView::fillVertexArray(ConnectionManager *remote, int length,
+void PView::fillVertexArray(ConnectionManager *remote, int length,
const char *bytes, int swap)
{
std::string name;
@@ -1402,7 +1402,7 @@ void PView::fillVertexArray(ConnectionManager *remote, int length,
if(!VertexArray::decodeHeader(length, bytes, swap, name, num, type, min, max,
numSteps, time, xmin, ymin, zmin, xmax, ymax, zmax))
return;
-
+
Msg::Debug("Filling vertex array (type %d) in view num %d", type, num);
SBoundingBox3d bbox(xmin, ymin, zmin, xmax, ymax, zmax);
@@ -1429,12 +1429,12 @@ void PView::fillVertexArray(ConnectionManager *remote, int length,
switch(type){
case 1:
- if(p->va_points) delete p->va_points;
+ if(p->va_points) delete p->va_points;
p->va_points = new VertexArray(1, 100);
p->va_points->fromChar(length, bytes, swap);
break;
- case 2:
- if(p->va_lines) delete p->va_lines;
+ case 2:
+ if(p->va_lines) delete p->va_lines;
p->va_lines = new VertexArray(2, 100);
p->va_lines->fromChar(length, bytes, swap);
break;
@@ -1453,7 +1453,7 @@ void PView::fillVertexArray(ConnectionManager *remote, int length,
p->va_ellipses = new VertexArray(4, 100);
p->va_ellipses->fromChar(length, bytes, swap);
break;
- default:
+ default:
Msg::Error("Cannot fill vertex array of type %d", type);
return;
}
--
GitLab