00001
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
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
00188 void debug_print( std::ostream &os ) const;
00189 };
00190
00191 int32_t _nodecount;
00192 int32_t _dof;
00193 Node2DoF _n2d;
00194 CRowMatrix *_fd_mat;
00195 Vector *_fd_vec;
00197 const Geometry *_g;
00198 mutable CRowMatrix *_fd_mat2;
00199 mutable Vector *_fd_vec2;
00200 mutable Vector *_fd_vec3;
00202 int32_t _neumann_order;
00203 bool _smooth_solid;
00205 plasma_mode_e _plasma;
00207 double _rhoe;
00208 double _Te;
00209 double _Up;
00211 std::vector<double> _rhoi;
00213 std::vector<double> _Ei;
00216 double _force_pot;
00218 bool (*_force_pot_func)(double,double,double);
00219 bool (*_init_plasma_func)(double,double,double);
00221 Solver *_solver;
00224 void set_link( CRowMatrix &A, Vector &B,
00225 int32_t a, int32_t b, double val );
00226
00227 void add_initial_plasma( int32_t i, int32_t j, int32_t k,
00228 CRowMatrix &A, Vector &B, Node2DoF &n2d );
00229
00230 void add_forced_pot( int32_t i, int32_t j, int32_t k,
00231 CRowMatrix &A, Vector &B, Node2DoF &n2d );
00232
00233 void add_vacuum_node( int32_t i, int32_t j, int32_t k,
00234 CRowMatrix &A, Vector &B, Node2DoF &n2d );
00235
00236 void add_neumann_node( signed char a, int32_t i, int32_t j, int32_t k,
00237 CRowMatrix &A, Vector &B, Node2DoF &n2d );
00238
00239 void add_solid_edge_node( signed char a, int32_t i, int32_t j, int32_t k,
00240 CRowMatrix &A, Vector &B, Node2DoF &n2d );
00241
00242 void clear_problem( void );
00243
00244 public:
00245
00246
00247
00248
00249
00252 EpotProblem();
00253
00256 EpotProblem( std::istream &s );
00257
00260 ~EpotProblem();
00261
00262
00263
00264
00265
00270 void set_neumann_order( int32_t order );
00271
00276 void enable_smooth_solids( bool enable );
00277
00285 void set_forced_potential_volume( double force_pot,
00286 bool (*force_pot_func)(double,double,double) );
00287
00292 void set_initial_plasma( double Up,
00293 bool (*plasma_func)(double,double,double) );
00294
00301 void set_pexp_plasma( double rhoe, double Te, double Up );
00302
00308 void set_nsimp_initial_plasma( bool (*plasma_func)(double,double,double) );
00309
00323 void set_nsimp_plasma( double rhop, double Ep,
00324 std::vector<double> rhoi, std::vector<double> Ei );
00325
00330 void construct( const Geometry &g );
00331
00334 void set_solver( Solver &s );
00335
00343 void solve( ScalarField &epot, const ScalarField &scharge ) const;
00344
00345
00346
00347
00348
00356 void get_vecmat( const Matrix **A, const Vector **B ) const;
00357
00365 void get_resjac( const Matrix **J, const Vector **R, const Vector &X ) const;
00366
00369 bool linear( void ) const;
00370
00371
00372
00373
00374
00377 int get_dof( void ) const { return( _dof ); }
00378
00381 void debug_print( std::ostream &os ) const;
00382
00385 void save( std::ostream &s ) const;
00386 };
00387
00388
00389 #endif
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405
00406
00407
00408
00409