LORENE
sol_elliptic_only_zec.C
1 /*
2  * Copyright (c) 2004 Philippe Grandclement
3  *
4  * This file is part of LORENE.
5  *
6  * LORENE is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License version 2
8  * as published by the Free Software Foundation.
9  *
10  * LORENE is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with LORENE; if not, write to the Free Software
17  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18  *
19  */
20 
21 char sol_elliptic_only_zec_C[] = "$Header: /cvsroot/Lorene/C++/Source/Non_class_members/PDE/sol_elliptic_only_zec.C,v 1.4 2014/10/13 08:53:30 j_novak Exp $" ;
22 
23 /*
24  * $Id: sol_elliptic_only_zec.C,v 1.4 2014/10/13 08:53:30 j_novak Exp $
25  * $Log: sol_elliptic_only_zec.C,v $
26  * Revision 1.4 2014/10/13 08:53:30 j_novak
27  * Lorene classes and functions now belong to the namespace Lorene.
28  *
29  * Revision 1.3 2014/10/06 15:16:10 j_novak
30  * Modified #include directives to use c++ syntax.
31  *
32  * Revision 1.2 2004/08/24 09:14:44 p_grandclement
33  * Addition of some new operators, like Poisson in 2d... It now requieres the
34  * GSL library to work.
35  *
36  * Also, the way a variable change is stored by a Param_elliptic is changed and
37  * no longer uses Change_var but rather 2 Scalars. The codes using that feature
38  * will requiere some modification. (It should concern only the ones about monopoles)
39  *
40  * Revision 1.1 2004/06/22 08:49:58 p_grandclement
41  * Addition of everything needed for using the logarithmic mapping
42  *
43  *
44  * $Header $
45  *
46  */
47 
48 // Header C :
49 #include <cstdlib>
50 #include <cmath>
51 
52 // Headers Lorene :
53 #include "tbl.h"
54 #include "mtbl_cf.h"
55 #include "map.h"
56 #include "param_elliptic.h"
57 
58 
59  //----------------------------------------------
60  // Version Mtbl_cf
61  //----------------------------------------------
62 
63 
64 
65 namespace Lorene {
66 Mtbl_cf elliptic_solver_only_zec (const Param_elliptic& ope_var, const Mtbl_cf& source, double val) {
67  // Verifications d'usage sur les zones
68  int nz = source.get_mg()->get_nzone() ;
69  assert (nz>1) ;
70  assert (source.get_mg()->get_type_r(nz-1) == UNSURR) ;
71 
72  // donnees sur la zone
73  int nr, nt, np ;
74 
75  //Rangement des valeurs intermediaires
76  Tbl *so ;
77  Tbl *sol_hom ;
78  Tbl *sol_part ;
79 
80 
81  // Rangement des solutions, avant raccordement
82  Mtbl_cf solution_part(source.get_mg(), source.base) ;
83  Mtbl_cf solution_hom_un(source.get_mg(), source.base) ;
84  Mtbl_cf solution_hom_deux(source.get_mg(), source.base) ;
85  Mtbl_cf resultat(source.get_mg(), source.base) ;
86 
87  solution_part.annule_hard() ;
88  solution_hom_un.annule_hard() ;
89  solution_hom_deux.annule_hard() ;
90  resultat.annule_hard() ;
91 
92  // Computation of the SP and SH's in the ZEC ...
93  int conte_start = 0 ;
94  for (int l=0 ; l<nz-1 ; l++) {
95  nt = source.get_mg()->get_nt(l) ;
96  np = source.get_mg()->get_np(l) ;
97  for (int k=0 ; k<np+1 ; k++)
98  for (int j=0 ; j<nt ; j++)
99  conte_start ++ ;
100  }
101  int conte = conte_start ;
102 
103  int zone = nz-1 ;
104 
105  nr = source.get_mg()->get_nr(zone) ;
106  nt = source.get_mg()->get_nt(zone) ;
107  np = source.get_mg()->get_np(zone) ;
108 
109  for (int k=0 ; k<np+1 ; k++)
110  for (int j=0 ; j<nt ; j++) {
111 
112  if (ope_var.operateurs[conte] != 0x0) {
113 
114  // Calcul de la SH
115  sol_hom = new Tbl(ope_var.operateurs[conte]->get_solh()) ;
116 
117  //Calcul de la SP
118  so = new Tbl(nr) ;
119  so->set_etat_qcq() ;
120  for (int i=0 ; i<nr ; i++)
121  so->set(i) = source(zone, k, j, i) ;
122 
123  sol_part = new Tbl(ope_var.operateurs[conte]->get_solp(*so)) ;
124 
125  // Rangement dans les tableaux globaux ;
126  for (int i=0 ; i<nr ; i++) {
127  solution_part.set(zone, k, j, i) = (*sol_part)(i) ;
128  if (sol_hom->get_ndim()==1)
129  solution_hom_un.set(zone, k, j, i) = (*sol_hom)(i) ;
130  else
131  {
132  solution_hom_un.set(zone, k, j, i) = (*sol_hom)(0,i) ;
133  solution_hom_deux.set(zone, k, j, i) = (*sol_hom)(1,i) ;
134  }
135  }
136 
137  delete so ;
138  delete sol_hom ;
139  delete sol_part ;
140 
141  }
142  conte ++ ;
143  }
144 
145 
146  //---------------------------------------------------------
147  // ON impose la bonne CL ... CASE ONLY SPHERICAL RIGHT NOW
148  //---------------------------------------------------------
149 
150  // C'est pas simple toute cette sombre affaire...
151  // Que le cas meme nombre de points dans chaque domaines...
152 
153  int start = conte_start ;
154  for (int k=0 ; k<1 ; k++)
155  for (int j=0 ; j<1 ; j++) {
156  if (ope_var.operateurs[start] != 0x0) {
157 
158 
159  // Valeurs en -1 :
160  double facteur = ((val-ope_var.F_minus(nz-1, k, j))/
161  ope_var.G_minus(nz-1) -
162  ope_var.operateurs[start]->val_sp_minus()) /
163  ope_var.operateurs[start]->val_sh_one_minus() ;
164 
165  // Zec
166  nr = source.get_mg()->get_nr(nz-1) ;
167  for (int i=0 ; i<nr ; i++)
168  resultat.set(nz-1,k,j,i) = solution_part(nz-1,k,j,i) +
169  facteur*solution_hom_un(nz-1,k,j,i) ;
170  }
171  start ++ ;
172  }
173 
174  return resultat;
175 }
176 
177 }
Lorene prototypes.
Definition: app_hor.h:64
const Mg3d * get_mg() const
Returns the Mg3d on which the Mtbl_cf is defined.
Definition: mtbl_cf.h:453
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448