diff --git a/common/wave/element.cpp b/common/wave/element.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..fda4158851d5b4cad66517290087865cfbb54cbb
--- /dev/null
+++ b/common/wave/element.cpp
@@ -0,0 +1,57 @@
+//GmshFWI Library
+#include "element.h"
+
+/*
+* WaveMultiField
+*/
+template<Physic T_Physic>
+void WaveMultiField<T_Physic>::initializeField(const WaveField<T_Physic>& field)
+{
+ _field = field;
+ for (auto it = _dofValues.begin(); it != _dofValues.end(); it++)
+ {
+ it->resize(_field.numberOfDofs());
+ }
+ _isInitialized = true;
+}
+
+template<Physic T_Physic>
+const WaveField<T_Physic>& WaveMultiField<T_Physic>::operator[](const unsigned int index) const
+{
+ _field.setAllVector(_dofValues[index], gmshfem::dofs::RawOrder::Hash);
+ return _field;
+}
+
+template<Physic T_Physic>
+void WaveMultiField<T_Physic>::setValues(const unsigned int index, gmshfem::algebra::Vector< std::complex< double > > &vector)
+{
+ _dofValues[index] = std::move(vector);
+}
+
+template<Physic T_Physic>
+void WaveMultiField<T_Physic>::write(std::string name) const
+{
+ for (unsigned int s = 0; s < size(); s++)
+ {
+ (*this)[s].write(name+"s"+std::to_string(s));
+ }
+}
+
+/*
+* WaveAuxilaryField
+*/
+template<Physic T_Physic>
+void WaveAuxilaryField<T_Physic>::assignValues(const std::vector< std::complex<double> > &values)
+{
+ unsigned int i = 0;
+ for(auto it = _wmf->_field.begin(); it != _wmf->_field.end(); it++)
+ {
+ if(it->first->type() == gmshfem::dofs::Type::Unknown)
+ {
+ _wmf->_dofValues[_index][i++] = values[it->first->numDof() - 1];
+ }
+ }
+}
+
+template class WaveMultiField<Physic::acoustic>;
+template class WaveAuxilaryField<Physic::acoustic>;
diff --git a/common/wave/element.h b/common/wave/element.h
index e69260daa5b271a2ff1bb71016ac42eaea033d65..1529b438439dd8058a63a616e8eb69b1213331ad 100644
--- a/common/wave/element.h
+++ b/common/wave/element.h
@@ -6,70 +6,34 @@
#include "../../specific/wave/element.h"
#include "../configuration.h"
+//Forward declaration
+template<Physic T_Physic> class WaveAuxilaryField;
+
/*
* WaveMultiField
*/
template<Physic T_Physic>
class WaveMultiField
{
-public:
+private:
mutable WaveField<T_Physic> _field;
bool _isInitialized;
std::vector< gmshfem::algebra::Vector< std::complex< double > > > _dofValues;
+public:
WaveMultiField(unsigned int n) : _field(), _isInitialized(false), _dofValues(n) {};
- void setField(const WaveField<T_Physic>& field)
- {
- _field = field;
- for (auto it = _dofValues.begin(); it != _dofValues.end(); it++)
- {
- it->resize(_field.numberOfDofs());
- }
- _isInitialized = true;
- };
-
+ void initializeField(const WaveField<T_Physic>& field);
bool isInitialized() const {return _isInitialized;};
- unsigned int size() const
- {
- return _dofValues.size();
- }
+ unsigned int size() const {return _dofValues.size();};
- const WaveField<T_Physic> &operator[](const unsigned int index) const
- {
- _field.setAllVector(_dofValues[index], gmshfem::dofs::RawOrder::Hash);
- return _field;
- }
+ const WaveField<T_Physic> &operator[](const unsigned int index) const;
- /*
- void setDofs(const std::vector< gmshfem::dofs::RawDof > &dofs)
- {
- _field->setAllDofVector(dofs, false);
- }
- */
+ void setValues(const unsigned int index, gmshfem::algebra::Vector< std::complex< double > > &vector);
- void setValues(const unsigned int index, gmshfem::algebra::Vector< std::complex< double > > &vector)
- {
- _dofValues[index] = std::move(vector);
- }
-/*
- std::vector< gmshfem::algebra::Vector< std::complex< double > > >::const_iterator begin() const
- {
- return _dofValues.begin();
- }
+ void write(std::string name) const;
- std::vector< gmshfem::algebra::Vector< std::complex< double > > >::const_iterator end() const
- {
- return _dofValues.end();
- }
-*/
- void write(std::string name) const
- {
- for (unsigned int s = 0; s < size(); s++)
- {
- (*this)[s].write(name+"s"+std::to_string(s));
- }
- };
+ friend class WaveAuxilaryField<T_Physic>;
};
/*
@@ -85,17 +49,8 @@ public:
WaveAuxilaryField(std::string name, gmshfem::domain::Domain domain, std::string gmodel, const wave::Discretization<T_Physic>& w_discret): WaveField<T_Physic>(name,domain,gmodel,w_discret), _wmf(nullptr), _index(0) {};
void setField(WaveMultiField<T_Physic>* wmf){_wmf=wmf;};
void setIndex(unsigned int index){_index=index;};
- virtual void assignValues(const std::vector< std::complex<double> > &values) override
- {
- unsigned int i = 0;
- for(auto it = _wmf->_field.begin(); it != _wmf->_field.end(); it++)
- {
- if(it->first->type() == gmshfem::dofs::Type::Unknown)
- {
- _wmf->_dofValues[_index][i++] = values[it->first->numDof() - 1];
- }
- }
- };
+ virtual void assignValues(const std::vector< std::complex<double> > &values) override;
};
+
#endif // H_COMMON_WAVE_ELEMENT
diff --git a/common/wave/equation/equation.cpp b/common/wave/equation/equation.cpp
index dd52426c0d171c2d86c09ea361728501c8f5cda4..94001d558ba83265470215814197078ef535ed29 100644
--- a/common/wave/equation/equation.cpp
+++ b/common/wave/equation/equation.cpp
@@ -104,7 +104,7 @@ DifferentialEquationInterface<T_Physic>::updateWithSystem(WaveMultiField<T_Physi
_formulation.removeTerms();
}
- if(!output.isInitialized()){output.setField(_v);};
+ if(!output.isInitialized()){output.initializeField(_v);};
_v.setField(&output);
for (unsigned int s = 0; s < _config->ns(); s++)
@@ -178,31 +178,33 @@ const WaveField<T_Physic>& DifferentialEquationInterface<T_Physic>::updateWithGr
template<Physic T_Physic>
void DifferentialEquationInterface<T_Physic>::updateGreen(const ModelStateEvaluator& m, const std::vector<unsigned int>& p_idx)
{
-// if(!_systemIsUpToDate)
-// {
-// _formulation.removeSystem();
-// _formulation.initSystem();
-//
-// _systemIsFactorized=false;
-// setLHS(m);
-// _formulation.pre();
-// _formulation.assemble();
-// _systemIsUpToDate=true;
-// _formulation.removeTerms();
-// }
-//
-// for (auto it = p_idx.begin(); it != p_idx.end(); it++)
-// {
-// setGreenRHS(*it);
-// _formulation.assemble();
-// _formulation.removeTerms();
-//
-// _formulation.solve(_systemIsFactorized);
-// _systemIsFactorized=true;
-// _formulation.setRHSToZero();
-// _g[*it] = _v;
-// //_v.getVector(_dofg[*it]);
-// }
+ if(!_systemIsUpToDate)
+ {
+ _formulation.removeSystem();
+ _formulation.initSystem();
+
+ _systemIsFactorized=false;
+ setLHS(m);
+ _formulation.pre();
+ _formulation.assemble();
+ _systemIsUpToDate=true;
+ _formulation.removeTerms();
+ }
+
+ if(!_g.isInitialized()){_g.initializeField(_v);};
+ _v.setField(&_g);
+
+ for (auto it = p_idx.begin(); it != p_idx.end(); it++)
+ {
+ _v.setIndex(*it);
+ setGreenRHS(*it);
+ _formulation.assemble();
+ _formulation.removeTerms();
+
+ _formulation.solve(_systemIsFactorized);
+ _systemIsFactorized=true;
+ _formulation.setRHSToZero();
+ }
}
template<Physic T_Physic>