IBSimu 1.0.4
|
00001 00005 /* Copyright (c) 2005-2010 Taneli Kalvas. All rights reserved. 00006 * 00007 * You can redistribute this software and/or modify it under the terms 00008 * of the GNU General Public License as published by the Free Software 00009 * Foundation; either version 2 of the License, or (at your option) 00010 * any later version. 00011 * 00012 * This library is distributed in the hope that it will be useful, but 00013 * WITHOUT ANY WARRANTY; without even the implied warranty of 00014 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00015 * General Public License for more details. 00016 * 00017 * You should have received a copy of the GNU General Public License 00018 * along with this library (file "COPYING" included in the package); 00019 * if not, write to the Free Software Foundation, Inc., 51 Franklin 00020 * Street, Fifth Floor, Boston, MA 02110-1301 USA 00021 * 00022 * If you have questions about your rights to use or distribute this 00023 * software, please contact Berkeley Lab's Technology Transfer 00024 * Department at TTD@lbl.gov. Other questions, comments and bug 00025 * reports should be sent directly to the author via email at 00026 * taneli.kalvas@jyu.fi. 00027 * 00028 * NOTICE. This software was developed under partial funding from the 00029 * U.S. Department of Energy. As such, the U.S. Government has been 00030 * granted for itself and others acting on its behalf a paid-up, 00031 * nonexclusive, irrevocable, worldwide license in the Software to 00032 * reproduce, prepare derivative works, and perform publicly and 00033 * display publicly. Beginning five (5) years after the date 00034 * permission to assert copyright is obtained from the U.S. Department 00035 * of Energy, and subject to any subsequent five (5) year renewals, 00036 * the U.S. Government is granted for itself and others acting on its 00037 * behalf a paid-up, nonexclusive, irrevocable, worldwide license in 00038 * the Software to reproduce, prepare derivative works, distribute 00039 * copies to the public, perform publicly and display publicly, and to 00040 * permit others to do so. 00041 */ 00042 00043 #ifndef EPOT_PROBLEM_HPP 00044 #define EPOT_PROBLEM_HPP 1 00045 00046 00047 #include <iostream> 00048 #include <stdint.h> 00049 #include "problem.hpp" 00050 #include "solver.hpp" 00051 #include "scalarfield.hpp" 00052 #include "geometry.hpp" 00053 #include "vec3d.hpp" 00054 00055 00068 enum plasma_mode_e {PLASMA_NONE = 0, PLASMA_PEXP_INITIAL, PLASMA_PEXP, 00069 PLASMA_NSIMP_INITIAL, PLASMA_NSIMP}; 00070 00071 #define PLASMA_INITIAL PLASMA_PEXP_INITIAL 00072 00142 class EpotProblem : public Problem { 00143 00153 class Node2DoF { 00154 Int3D _size; 00155 int32_t *_n2d; 00157 public: 00158 00159 Node2DoF() : _size(0), _n2d(0) {} 00160 Node2DoF( Int3D size ) : _size(size) { 00161 _n2d = new int32_t[_size[0]*_size[1]*_size[2]]; 00162 } 00163 ~Node2DoF() { delete _n2d; } 00164 00165 void resize( Int3D size ) { 00166 _size = size; 00167 if( _n2d ) 00168 delete _n2d; 00169 _n2d = new int32_t[_size[0]*_size[1]*_size[2]]; 00170 } 00171 00172 int32_t &operator()( int i ) 00173 { return( _n2d[i] ); } 00174 int32_t &operator()( int i, int j ) 00175 { return( _n2d[i+j*_size[0]] ); } 00176 int32_t &operator()( int i, int j, int k ) 00177 { return( _n2d[i+j*_size[0]+k*_size[0]*_size[1]] ); } 00178 00179 const int32_t &operator()( int i ) const 00180 { return( _n2d[i] ); } 00181 const int32_t &operator()( int i, int j ) const 00182 { return( _n2d[i+j*_size[0]] ); } 00183 const int32_t &operator()( int i, int j, int k ) const 00184 { return( _n2d[i+j*_size[0]+k*_size[0]*_size[1]] ); } 00185 00186 void debug_print( void ) const; 00187 }; 00188 00189 int32_t _nodecount; 00190 int32_t _dof; 00191 Node2DoF _n2d; 00192 CRowMatrix *_fd_mat; 00193 Vector *_fd_vec; 00195 const Geometry *_g; 00196 mutable CRowMatrix *_fd_mat2; 00197 mutable Vector *_fd_vec2; 00198 mutable Vector *_fd_vec3; 00200 int32_t _neumann_order; 00201 bool _smooth_solid; 00203 plasma_mode_e _plasma; 00205 double _rhoe; 00206 double _Te; 00207 double _Up; 00209 std::vector<double> _rhoi; 00211 std::vector<double> _Ei; 00214 double _force_pot; 00216 bool (*_force_pot_func)(double,double,double); 00217 bool (*_init_plasma_func)(double,double,double); 00219 Solver *_solver; 00222 void set_link( CRowMatrix &A, Vector &B, 00223 int32_t a, int32_t b, double val ); 00224 00225 void add_initial_plasma( int32_t i, int32_t j, int32_t k, 00226 CRowMatrix &A, Vector &B, Node2DoF &n2d ); 00227 00228 void add_forced_pot( int32_t i, int32_t j, int32_t k, 00229 CRowMatrix &A, Vector &B, Node2DoF &n2d ); 00230 00231 void add_vacuum_node( int32_t i, int32_t j, int32_t k, 00232 CRowMatrix &A, Vector &B, Node2DoF &n2d ); 00233 00234 void add_neumann_node( signed char a, int32_t i, int32_t j, int32_t k, 00235 CRowMatrix &A, Vector &B, Node2DoF &n2d ); 00236 00237 void add_solid_edge_node( signed char a, int32_t i, int32_t j, int32_t k, 00238 CRowMatrix &A, Vector &B, Node2DoF &n2d ); 00239 00240 void clear_problem( void ); 00241 00242 public: 00243 00244 /* ************************************** * 00245 * Constructors and destructor * 00246 * ************************************** */ 00247 00250 EpotProblem(); 00251 00254 EpotProblem( std::istream &s ); 00255 00258 ~EpotProblem(); 00259 00260 /* ************************************** * 00261 * Problem constructing and solving * 00262 * ************************************** */ 00263 00268 void set_neumann_order( int32_t order ); 00269 00274 void enable_smooth_solids( bool enable ); 00275 00283 void set_forced_potential_volume( double force_pot, 00284 bool (*force_pot_func)(double,double,double) ); 00285 00290 void set_initial_plasma( double Up, 00291 bool (*plasma_func)(double,double,double) ); 00292 00295 void set_pexp_plasma( double rhoe, double Te, double Up ); 00296 00302 void set_nsimp_initial_plasma( bool (*plasma_func)(double,double,double) ); 00303 00317 void set_nsimp_plasma( double rhop, double Ep, 00318 std::vector<double> rhoi, std::vector<double> Ei ); 00319 00324 void construct( const Geometry &g ); 00325 00328 void set_solver( Solver &s ); 00329 00337 void solve( ScalarField &epot, const ScalarField &scharge ) const; 00338 00339 /* ************************************** * 00340 * Solver interface * 00341 * ************************************** */ 00342 00350 void get_vecmat( const Matrix **A, const Vector **B ) const; 00351 00359 void get_resjac( const Matrix **J, const Vector **R, const Vector &X ) const; 00360 00363 bool linear( void ) const; 00364 00365 /* ************************************** * 00366 * Misc * 00367 * ************************************** */ 00368 00371 int get_dof( void ) const { return( _dof ); } 00372 00375 void debug_print( void ) const; 00376 00379 void save( std::ostream &s ) const; 00380 }; 00381 00382 00383 #endif 00384 00385 00386 00387 00388 00389 00390 00391 00392 00393 00394 00395 00396 00397 00398 00399 00400 00401 00402