20 #include "system_of_eqs.hpp"
22 #include "tensor_impl.hpp"
24 #include "exceptions.hpp"
34 for (
int i=0 ; i<
neq_int ; i++) {
35 errors.
set(i) = fabs(sec(pos)) ;
39 for (
int i=0 ; i<
neq ; i++) {
41 for (
int j=0 ; j<
eq[i]->get_n_cond_tot() ; j++) {
42 if (fabs(sec(pos)) > max)
43 max = fabs(sec(pos)) ;
57 for (
int i=0 ; i<
ndef ; i++)
58 def[i]->compute_res() ;
61 for (
int i=0 ; i<
neq ; i++)
69 for (
int i=0 ; i<
neq ; i++)
84 for (
int i=0 ; i<
neq_int ; i++) {
89 for (
int i=0 ; i<
neq ; i++)
90 eq[i]->export_val(conte,
results.set_data(), res, pos_res) ;
101 for (
int i=0 ; i<
ndef ; i++)
102 def[i]->compute_res() ;
105 for (
int i=0 ; i<
neq ; i++)
106 eq[i]->apply(conte,
results.set_data()) ;
109 cerr <<
"Number of conditions unknown ; call sec_member first" << endl ;
117 for (
int i=0 ; i<
neq_int ; i++) {
118 res.
set(pos_res) =
eq_int[i]->get_der() ;
122 for (
int i=0 ; i<
neq ; i++)
123 eq[i]->export_der(conte,
results.set_data(), res, pos_res) ;
134 bool is_var_double = false ;
143 cst[i]->set_der_zero() ;
144 for (
int i=0 ; i<
nterm ; i++)
145 term[i]->set_der_zero() ;
153 is_var_double = true ;
163 for (
int i=0 ; i<
nterm ; i++) {
164 int dom =
term[i]->get_dom() ;
170 catch(Unknown_base_error & e) {
171 std::cerr <<
"Error in System_of_eqs[rank=" <<
mpi_proc_rank <<
"]::do_col_J(cc = " << cc
175 std::cerr <<
" with auxi=term[i]->get_val_t(), where i="<< i ;
176 std::cerr <<
", dom=" << dom <<
" and conte=" << conte << std::endl;
181 if ((!auxi(auxi.
indices(j))(dom).check_if_zero()) && (!auxi(auxi.
indices(j))(dom).get_base().is_def()))
183 if ((zedom==-1) && (!auxi(auxi.
indices(j))(dom).check_if_zero()))
187 term[i]->set_der_t(auxi +
term[i]->get_der_t()) ;
196 for (
int i=0 ; i<
ndef ; i++)
197 def[i]->compute_res() ;
200 for (
int i=0 ; i<
neq ; i++) {
201 if ((is_var_double) || (
eq[i]->take_into_account(zedom)) || (
eq[i]->take_into_account(zedoms(0))) || (
eq[i]->take_into_account(zedoms(1)))) {
202 eq[i]->apply(conte,
results.set_data()) ;
205 conte +=
eq[i]->n_ope ;
209 cerr <<
"Number of conditions unknown ; call sec_member first" << endl ;
217 for (
int i=0 ; i<
neq_int ; i++) {
218 res.
set(pos_res) =
eq_int[i]->get_der() ;
222 for (
int i=0 ; i<
neq ; i++) {
223 if ((is_var_double) || (
eq[i]->take_into_account(zedom))|| (
eq[i]->take_into_account(zedoms(0))) || (
eq[i]->take_into_account(zedoms(1)))) {
224 eq[i]->export_der(conte,
results.set_data(), res, pos_res) ;
228 pos_res +=
eq[i]->n_cond_tot ;
229 conte +=
eq[i]->n_ope ;
238 for (
int i=0 ; i<
nterm ; i++) {
239 int dom =
term[i]->get_dom() ;
241 for (
int d=0 ; d<zedoms.
get_size(0) ; d++)
247 term[i]->set_der_zero() ;
251 if (
cst[i]->get_type_data() == TERM_T) {
252 int dom =
cst[i]->get_dom() ;
254 for (
int d=0 ; d<zedoms.
get_size(0) ; d++)
260 cst[i]->set_der_zero() ;
reference set(const Index &pos)
Read/write of an element.
int get_size(int i) const
Returns the size of a given dimension.
virtual void update_term_eq(Term_eq *so) const
Update the value of a field, after the shape of the Domain has been changed by the system.
virtual void affecte_tau_one_coef(Tensor &so, int dom, int cc, int &pos_cf) const
Sets at most one coefficient of a Tensor to 1.
virtual void update()
Updates the derived quantities (Christoffels etc..) This is done only for the ones that are needed,...
Val_domain & set_domain(int)
Read/write of a particular Val_domain.
const Domain * get_domain(int i) const
returns a pointer on the domain.
virtual void affecte_coef_to_variable_domains(int &conte, int cc, Array< int > &doms) const
The variation of the functions describing the shape of the Domain are affected from the unknowns of t...
MMPtr_array< Ope_def > def
Pointers on the definition (i.e. on the Ope_def that is needed to compute the result).
Array< double > sec_member()
Computes the second member of the Newton-Raphson equations.
int nvar_double
Number of unknowns that are numbers (i.e. not fields)
MMPtr_array< Term_eq > results
Pointers on the residual of the various equations.
int nterm
Number of Term_eq corresponding to the unknown fields.
Array< double > do_col_J(int i)
Computes one column of the Jacobian.
int mpi_proc_rank
Sylvain's stuff.
MMPtr_array< Term_eq > term_double
Pointers on the Term_eq corresponding to the unknowns that are numbers.
int nbr_unknowns
Number of unknowns (basically the number of coefficients of all the unknown fields,...
int ndom
Number of domains used.
int neq_int
Number of integral equations (i.e. which are doubles)
void update_terms_from_variable_domains(const Array< int > &zedoms)
Updates the variations of the Term_eq that comes from the fact that some Domains are variable (i....
MMPtr_array< Equation > eq
Pointers onto the field equations.
virtual void compute_nbr_of_conditions()
Sylvain's stuff.
MMPtr_array< Term_eq > cst
Pointers on the Term_eq coming from the constants passed by the user.
Array< double > check_equations()
Computes the residual of all the equations.
void vars_to_terms()
Copies the various unknowns (doubles and Tensors) into their Term_eq counterparts.
int ndef
Number of definitions.
int nterm_cst
Number of Term_eq coming from the constants passed by the user.
int nbr_conditions
Total number of conditions (the number of coefficients of all the equations, once regularities are ta...
MMPtr_array< Term_eq > term
Pointers on the Term_eq corresponding to the unknown fields.
int neq
Number of field equations.
int dom_max
Highest domain number.
const Space & espace
Associated Space.
Array< double > do_JX(const Array< double > &xx)
Computes the product where is the Jacobian of the system.
MMPtr_array< Eq_int > eq_int
Pointers onto the integral equations.
void xx_to_ders(const Array< double > &vder)
Sets the values the variation of the fields.
int dom_min
Smallest domain number.
Metric * met
Pointer on the associated Metric, if defined.
Scalar & set(const Array< int > &ind)
Returns the value of a component (read/write version).
int get_n_comp() const
Returns the number of stored components.
virtual Array< int > indices(int pos) const
Gives the values of the indices corresponding to a location in the array used for storage of the comp...
Base_spectral & set_base()
Sets the basis of decomposition.