particledatabase.hpp

Go to the documentation of this file.
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 PARTICLEDATABASE_HPP
00044 #define PARTICLEDATABASE_HPP 1
00045 
00046 
00047 #include <vector>
00048 #include "timer.hpp"
00049 #include "ibsimu.hpp"
00050 #include "trajectory.hpp"
00051 #include "particles.hpp"
00052 #include "particleiterator.hpp"
00053 #include "trajectorydiagnostics.hpp"
00054 
00055 
00056 
00057 /* ******************************************************************************************* *
00058  * ParticleDataBase classes                                                                    *
00059  * ******************************************************************************************* */
00060 
00061 
00069 class ParticleDataBase {
00070 
00071 protected:
00072 
00073     uint32_t       _threadcount; 
00074     double         _epsabs;      
00075     double         _epsrel;      
00076     bool           _polyint;     
00077     uint32_t       _maxsteps;    
00078     double         _maxt;        
00079     uint32_t       _trajdiv;     
00081     bool           _mirror[6];   
00083     double         _rhosum;      
00085     uint32_t       _end_time;    
00086     uint32_t       _end_step;    
00087     uint32_t       _end_out;     
00089     uint32_t       _end_coll;    
00091     uint32_t       _end_baddef;  
00092     uint32_t       _sum_steps;   
00094     int            _iteration;   
00096     bool           _nsimp;
00097     const ScalarField *_epot;
00098     double         _phi_plasma;
00099 
00102     ParticleDataBase()
00103         : _threadcount(1), _epsabs(1e-6), _epsrel(1e-6), _polyint(true), _maxsteps(1000), 
00104           _maxt(1e-3), _trajdiv(1), _rhosum(0.0), _end_time(0), _end_step(0), _end_out(0), 
00105           _end_coll(0), _end_baddef(0), _sum_steps(0), _iteration(-1), _nsimp(false), 
00106           _epot(NULL), _phi_plasma(0.0) {
00107         _mirror[0] = false;
00108         _mirror[1] = false;
00109         _mirror[2] = false;
00110         _mirror[3] = false;
00111         _mirror[4] = false;
00112         _mirror[5] = false;
00113     }
00114 
00115 public:
00116 
00117 /* ************************************** *
00118  * Constructors and destructor            *
00119  * ************************************** */
00120 
00123     virtual ~ParticleDataBase() {}
00124 
00125 /* ****************************************** *
00126  * Particle iteration settings and statictics *
00127  * ****************************************** */
00128 
00131     void set_thread_count( uint32_t threadcount ) {
00132         if( threadcount <= 0 )
00133             throw( Error( ERROR_LOCATION, "invalid parameter" ) );
00134         _threadcount = threadcount;
00135     }
00136 
00142     void set_accuracy( double epsabs, double epsrel ) {
00143         _epsabs = epsabs;
00144         _epsrel = epsrel;
00145     }
00146 
00147     void enable_nsimp_plasma_threshold( const ScalarField *epot, double phi_plasma ) {
00148         _nsimp = true;
00149         _epot = epot;
00150         _phi_plasma = phi_plasma;
00151     }
00152 
00157     void set_polyint( bool polyint ) {
00158         _polyint = polyint;
00159     }
00160     
00166     bool get_polyint( void ) const {
00167         return( _polyint );
00168     }
00169     
00174     void set_max_steps( uint32_t maxsteps ) {
00175         if( maxsteps <= 0 )
00176             throw( Error( ERROR_LOCATION, "invalid parameter" ) );
00177         _maxsteps = maxsteps;
00178     }
00179 
00184     void set_max_time( double maxt ) {
00185         if( maxt <= 0.0 )
00186             throw( Error( ERROR_LOCATION, "invalid parameter" ) );
00187         _maxt = maxt;
00188     }
00189 
00196     void set_save_trajectories( uint32_t div ) {
00197         _trajdiv = div;
00198     }
00199 
00204     void set_mirror( const bool mirror[6] ) {
00205         _mirror[0] = mirror[0];
00206         _mirror[1] = mirror[1];
00207         _mirror[2] = mirror[2];
00208         _mirror[3] = mirror[3];
00209         _mirror[4] = mirror[4];
00210         _mirror[5] = mirror[5];
00211     }
00212 
00217     void get_mirror( bool mirror[6] ) const {
00218         mirror[0] = _mirror[0];
00219         mirror[1] = _mirror[1];
00220         mirror[2] = _mirror[2];
00221         mirror[3] = _mirror[3];
00222         mirror[4] = _mirror[4];
00223         mirror[5] = _mirror[5];
00224     }
00225 
00226     int get_iteration_number( void ) const {
00227         return( _iteration );
00228     }
00229 
00243     double get_rhosum( void ) {
00244         return( _rhosum );
00245     }
00246 
00247 /* ************************************** *
00248  * Information and queries                *
00249  * ************************************** */
00250 
00253     virtual size_t size( void ) const = 0;
00254 
00257     virtual ParticleBase &particle( uint32_t i ) = 0;
00258 
00261     virtual const ParticleBase &particle( uint32_t i ) const = 0;
00262 
00265     virtual size_t traj_size( uint32_t i ) const = 0;
00266     
00269     virtual void trajectory_point( double &t, Vec3D &loc, Vec3D &vel, uint32_t i, uint32_t j ) const = 0;
00270 
00274     virtual void trajectories_at_plane( TrajectoryDiagnosticData &tdata, 
00275                                         coordinate_axis_e axis,
00276                                         double val,
00277                                         const std::vector<trajectory_diagnostic_e> &diagnostics ) const = 0;
00278 
00279 /* ************************************** *
00280  * Particle and trajectory clearing       *
00281  * ************************************** */
00282 
00288     virtual void clear( void ) = 0;
00289 
00295     virtual void clear_trajectories( void ) = 0;
00296 
00297 /* ************************************** *
00298  * Particle definition                    *
00299  * ************************************** */
00300 
00303     virtual void reserve( size_t size ) = 0;
00304 
00305 /* ************************************** *
00306  * Debugging, plotting and saving         *
00307  * ************************************** */
00308 
00309     virtual void debug_print( void ) const = 0;
00310 
00311 };
00312 
00313 
00328 template<class PP> class ParticleDataBasePP : public ParticleDataBase {
00329 
00332     static void add_diagnostics( TrajectoryDiagnosticData &tdata, const PP &x, 
00333                                  const Particle<PP> &p, int crd ) {
00334         //std::cout << "add_diagnostics():\n";
00335         for( size_t a = 0; a < tdata.diag_size(); a++ ) {
00336             //std::cout << "  diagnostic[" << a << "] = " << tdata.diagnostic(a) << "\n";
00337             
00338             double data = 0.0;
00339             switch( tdata.diagnostic( a ) ) {
00340             case DIAG_NONE:
00341                 data = 0.0;
00342                 break;
00343             case DIAG_T:
00344                 data = x[0];
00345                 break;
00346             case DIAG_X:
00347                 data = x[1];
00348                 break;
00349             case DIAG_VX:
00350                 data = x[2];
00351                 break;
00352             case DIAG_Y:
00353             case DIAG_R:
00354                 data = x[3];
00355                 break;
00356             case DIAG_VY:
00357             case DIAG_VR:
00358                 data = x[4];
00359                 break;
00360             case DIAG_Z:
00361                 data = x[5];
00362                 break;
00363             case DIAG_VZ:
00364                 data = x[6];
00365                 break;
00366             case DIAG_W:
00367                 data = x[5];
00368                 break;
00369             case DIAG_VTHETA:
00370                 data = x[5]*x[3];
00371                 break;
00372             case DIAG_XP:
00373                 data = x[2]/x[2*crd+2];
00374                 break;
00375             case DIAG_YP:
00376             case DIAG_RP:
00377                 data = x[4]/x[2*crd+2];
00378                 break;
00379             case DIAG_AP:
00380                 data = x[3]*x[5]/x[2*crd+2];
00381                 break;
00382             case DIAG_ZP:
00383                 data = x[6]/x[2*crd+2];
00384                 break;
00385             case DIAG_CURR:
00386                 data = p.IQ();
00387                 break;
00388             case DIAG_QM:
00389                 data = p.qm();
00390                 break;
00391             case DIAG_EK:
00392                 Vec3D velocity = x.velocity();
00393                 data = velocity.norm2();
00394                 break;
00395             }
00396             //std::cout << "  adding data = " << data << "\n";
00397             tdata.add_data( a, data );
00398         }
00399     }
00400 
00401 protected:
00402 
00403     std::vector<Particle<PP> > _particles;      
00407     ParticleDataBasePP() {}
00408 
00409 public:
00410 
00411 /* ************************************** *
00412  * Constructors and destructor            *
00413  * ************************************** */
00414 
00417     ~ParticleDataBasePP() {}
00418 
00419 /* ************************************** *
00420  * Information and queries                *
00421  * ************************************** */
00422 
00425     virtual size_t size( void ) const { return( _particles.size() ); }
00426 
00429     virtual Particle<PP> &particle( uint32_t i ) { return( _particles[i] ); }
00430 
00433     virtual const Particle<PP> &particle( uint32_t i ) const { return( _particles[i] ); }
00434     
00437     virtual size_t traj_size( uint32_t i ) const { return( _particles[i].traj_size() ); }
00438     
00441     virtual void trajectory_point( double &t, Vec3D &loc, Vec3D &vel, uint32_t i, uint32_t j ) const {
00442         PP x = _particles[i].traj(j); 
00443         t = x[0];
00444         loc = x.location();
00445         vel = x.velocity();
00446     }
00447 
00451     virtual void trajectories_at_plane( TrajectoryDiagnosticData &tdata, 
00452                                         coordinate_axis_e axis,
00453                                         double val,
00454                                         const std::vector<trajectory_diagnostic_e> &diagnostics ) const {
00455 
00456         if( ibsimu.get_verbose_output() )
00457             std::cout << "Making trajectory diagnostics at " 
00458                       << coordinate_axis_string[axis] << " = " << val << "\n";
00459 
00460         // Check query
00461         switch( PP::geom_mode() ) {
00462         case MODE_1D:
00463             throw( Error( ERROR_LOCATION, "unsupported dimension number" ) );
00464             break;
00465         case MODE_2D:
00466             if( axis == AXIS_R || axis == AXIS_Z )
00467                 throw( Error( ERROR_LOCATION, "nonexistent axis" ) );
00468             break;
00469         case MODE_CYL:
00470             if( axis == AXIS_Y || axis == AXIS_Z )
00471                 throw( Error( ERROR_LOCATION, "nonexistent axis" ) );
00472             break;
00473         case MODE_3D:
00474             if( axis == AXIS_R )
00475                 throw( Error( ERROR_LOCATION, "nonexistent axis" ) );
00476             break;
00477         default:
00478             throw( Error( ERROR_LOCATION, "unsupported dimension number" ) );
00479         }
00480 
00481         // Check diagnostics query validity
00482         for( size_t a = 0; a < diagnostics.size(); a++ ) {
00483             if( diagnostics[a] == DIAG_NONE )
00484                 throw( Error( ERROR_LOCATION, "invalid diagnostics query \'DIAG_NONE\'" ) );
00485             else if( PP::geom_mode() != MODE_CYL && (diagnostics[a] == DIAG_R ||
00486                                                      diagnostics[a] == DIAG_VR ||
00487                                                      diagnostics[a] == DIAG_RP ||
00488                                                      diagnostics[a] == DIAG_W ||
00489                                                      diagnostics[a] == DIAG_VTHETA ||
00490                                                      diagnostics[a] == DIAG_AP) )
00491                 throw( Error( ERROR_LOCATION, "invalid diagnostics query for geometry type" ) );
00492         }
00493 
00494         // Prepare output vector
00495         tdata.clear();
00496         for( size_t a = 0; a < diagnostics.size(); a++ ) {
00497             tdata.add_data_column( diagnostics[a] );
00498         }
00499         
00500         // Set coordinate index
00501         int crd;
00502         switch( axis ) {
00503         case AXIS_X:
00504             crd = 0;
00505             break;
00506         case AXIS_Y:
00507         case AXIS_R:
00508             crd = 1;
00509             break;
00510         case AXIS_Z:
00511             crd = 2;
00512             break;
00513         default:
00514             throw( Error( ERROR_LOCATION, "unsupported axis" ) );
00515         }
00516 
00517         // Scan through particle trajectory points
00518         double Isum = 0.0;
00519         std::vector<PP> intsc;
00520         for( size_t a = 0; a < _particles.size(); a++ ) {
00521             size_t N = _particles[a].traj_size();
00522             if( N < 2 )
00523                 continue;
00524             PP x1 = _particles[a].traj(0);
00525             size_t nintsum = 0;
00526             for( size_t b = 1; b < N; b++ ) {
00527                 PP x2 = _particles[a].traj(b);
00528                 intsc.clear();
00529                 size_t nintsc = PP::trajectory_intersections_at_plane( intsc, crd, val, x1, x2 );
00530                 nintsum += nintsc;
00531                 for( size_t c = 0; c < nintsc; c++ ) {
00532                     Isum += _particles[a].IQ();
00533                     add_diagnostics( tdata, intsc[c], _particles[a], crd );
00534                 }
00535 
00536                 x1 = x2;
00537             }
00538         }
00539 
00540         if( ibsimu.get_verbose_output() ) {
00541             std::cout << "  number of trajectories = " << tdata.traj_size() << "\n";
00542             if( PP::geom_mode() == MODE_2D )
00543                 std::cout << "  total current = " << Isum << " A/m\n";
00544             else
00545                 std::cout << "  total current = " << Isum << " A\n";
00546         }
00547     }
00548 
00549 /* ************************************** *
00550  * Particle and trajectory clearing       *
00551  * ************************************** */
00552 
00555     virtual void clear( void ) { 
00556         _particles.clear();
00557         _rhosum = 0.0;
00558     }
00559 
00565     virtual void clear_trajectories( void ) {
00566         for( uint32_t a = 0; a < _particles.size(); a++ )
00567             _particles[a].clear_trajectory();
00568     }
00569 
00570 /* ************************************** *
00571  * Particle definition                    *
00572  * ************************************** */
00573 
00576     virtual void reserve( size_t size ) { _particles.reserve( size ); }
00577 
00589     void add_particle( double IQ, double q, double m, const PP &x ) {
00590         _particles.push_back( Particle<PP>( IQ, CHARGE_E*q, MASS_U*m, x ) );
00591     }
00592 
00597     void add_particle( const Particle<PP> &pp ) {
00598         _particles.push_back( pp );
00599     }
00600 
00601 /* ************************************** *
00602  * Particle iterators                     *
00603  * ************************************** */
00604 
00612     virtual void iterate_trajectories( ScalarField &scharge, const Efield &efield, 
00613                                        const VectorField &bfield, const Geometry &g ) {
00614 
00615         ScalarField                         *schmap[_threadcount];
00616         std::vector<ParticleIterator<PP> *>  iterators;
00617 
00618         Timer t;
00619         if( ibsimu.get_verbose_output() )
00620             std::cout << "Calculating particle trajectories\n";
00621         _iteration++;
00622 
00623         // Check geometry mode
00624         if( g.geom_mode() != PP::geom_mode() )
00625             throw( Error( ERROR_LOCATION, "Differing geometry modes" ) );
00626 
00627         // Clear space charge
00628         scharge.clear();
00629 
00630         // Clear statistics
00631         _end_time   = 0;
00632         _end_step   = 0;
00633         _end_out    = 0;
00634         _end_coll   = 0;
00635         _end_baddef = 0;
00636         _sum_steps  = 0;
00637 
00638         // Check number of particles
00639         if( _particles.size() == 0 ) {
00640             std::cout << "  no particles to calculate\n";
00641             return;
00642         }
00643 
00644         // Make separate space charge maps for all threads and build iterators
00645         for( uint32_t a = 0; a < _threadcount; a++ ) {
00646             if( a == 0 ) schmap[a] = &scharge;
00647             else schmap[a] = new ScalarField( g );
00648             iterators.push_back( new ParticleIterator<PP>( PARTICLE_ITERATOR_ADAPTIVE, _epsabs, _epsrel, 
00649                                                            _polyint, _maxsteps, 
00650                                                            _maxt, _trajdiv, _mirror, schmap[a], 
00651                                                            &efield, &bfield, &g, &_particles[0] ) );
00652             if( _nsimp )
00653                 iterators[a]->enable_nsimp_plasma_threshold( _epot, _phi_plasma );
00654         }
00655 
00656         // Make Scheduler
00657         Scheduler<ParticleIterator<PP>,Particle<PP>,Error> scheduler( iterators );
00658 
00659         // Add problems
00660         for( size_t a = 0; a < _particles.size(); a++ )
00661             scheduler.add_problem( &_particles[a] );
00662 
00663         // Wait for completition
00664         scheduler.run();
00665         scheduler.finish();
00666 
00667         if( scheduler.is_error() ) {
00668             // Throw the error
00669             std::vector<Error> err;
00670             std::vector<Particle<PP> *> part;
00671             scheduler.get_errors( err, part );
00672             throw( err[0] );
00673         }
00674 
00675         // Combine separate space charge maps and collect
00676         // statistics. Free all allocated memory.
00677         for( uint32_t a = 0; a < _threadcount; a++ ) {
00678             if( a != 0 ) {
00679                 scharge += *schmap[a];
00680                 delete schmap[a];
00681             }
00682             uint32_t end_time, end_step, end_out, end_coll, end_baddef, sum_steps;
00683             iterators[a]->get_statistics( end_time, end_step, end_out, 
00684                                           end_coll, end_baddef, sum_steps );
00685             _end_time   += end_time;
00686             _end_step   += end_step;
00687             _end_out    += end_out;
00688             _end_coll   += end_coll;
00689             _end_baddef += end_baddef;
00690             _sum_steps  += sum_steps;
00691             delete iterators[a];
00692         }
00693 
00694         scharge_finalize( scharge );
00695         
00696         t.stop();
00697         if( ibsimu.get_verbose_output() ) {
00698             std::cout << "  Particle histories (" << _particles.size() << " total):\n";
00699             std::cout << "    time limited = " << _end_time << "\n";
00700             std::cout << "    step count limited = " << _end_step << "\n";
00701             std::cout << "    out of geometry = " << _end_out << "\n";
00702             std::cout << "    collided = " << _end_coll << "\n";
00703             std::cout << "    bad definitions = " << _end_baddef << "\n";
00704             std::cout << "    total steps = " << _sum_steps << "\n";
00705             std::cout << "    steps per particle (ave) = " << _sum_steps/(double)_particles.size() << "\n";
00706             std::cout << "  time used = " << t << "\n";
00707         }
00708     }
00709 
00716     virtual void step_particles( const Efield &efield, const Geometry &g, double dt ) {
00717 
00718     }
00719 
00720 /* ************************************** *
00721  * Debugging, plotting and saving         *
00722  * ************************************** */
00723 
00726     virtual void debug_print( void ) const {
00727         std::cout << "threadcount = " << _threadcount << "\n";
00728         std::cout << "epsabs = "      << _epsabs << "\n";
00729         std::cout << "epsrel = "      << _epsrel << "\n";
00730         std::cout << "maxsteps = "    << _maxsteps << "\n";
00731         std::cout << "maxt = "        << _maxt << "\n";
00732         std::cout << "trajdiv = "     << _trajdiv << "\n";
00733         std::cout << "mirror = (";
00734         for( uint32_t a = 0; a < 5; a++ )
00735             std::cout << _mirror[a] << ", ";
00736         std::cout << _mirror[5] << ")\n";
00737         
00738         for( uint32_t a = 0; a < _particles.size(); a++ ) {
00739             std::cout << "Particle " << a << ":\n";
00740             _particles[a].debug_print();
00741         }
00742     }
00743 
00744 };
00745 
00746 
00747 
00759 class ParticleDataBase2D : public ParticleDataBasePP<ParticleP2D> {
00760 
00761 
00762 public:
00763 
00764 /* ************************************** *
00765  * Constructors and destructor            *
00766  * ************************************** */
00767 
00770     ParticleDataBase2D() {}
00771 
00774     ~ParticleDataBase2D() {}
00775 
00776 /* ************************************** *
00777  * Particle beam definition               *
00778  * ************************************** */
00779 
00799     void add_2d_beam_with_energy( uint32_t N, double J, double q, double m, 
00800                                   double E, double Tp, double Tt, 
00801                                   double x1, double y1, double x2, double y2 );
00802 
00818     void add_2d_beam_with_velocity( uint32_t N, double J, double q, double m, 
00819                                     double v, double dvp, double dvt, 
00820                                     double x1, double y1, double x2, double y2 );
00821 
00838     void add_2d_KV_beam_with_emittance( uint32_t N, double I, double q, double m,
00839                                         double a, double b, double e,
00840                                         double Ex, double x0, double y0 );
00841 
00858     void add_2d_gaussian_beam_with_emittance( uint32_t N, double I, double q, double m,
00859                                               double a, double b, double e,
00860                                               double Ex, double x0, double y0 );
00861 };
00862 
00863 
00864 
00865 
00877 class ParticleDataBaseCyl : public ParticleDataBasePP<ParticlePCyl> {
00878 
00879 
00880 public:
00881 
00882 /* ************************************** *
00883  * Constructors and destructor            *
00884  * ************************************** */
00885 
00888     ParticleDataBaseCyl() {}
00889 
00892     ~ParticleDataBaseCyl () {}
00893 
00894 /* ************************************** *
00895  * Particle beam definition               *
00896  * ************************************** */
00897 
00917     void add_2d_beam_with_energy( uint32_t N, double J, double q, double m, 
00918                                   double E, double Tp, double Tt, 
00919                                   double x1, double y1, double x2, double y2 );
00920 
00936     void add_2d_beam_with_velocity( uint32_t N, double J, double q, double m, 
00937                                     double v, double dvp, double dvt, 
00938                                     double x1, double y1, double x2, double y2 );
00939 
00942     void add_2d_full_gaussian_beam( uint32_t N, double I, double q, double m,
00943                                     double Ex, double Tp, double Tt, 
00944                                     double x0, double dr );
00945 
00964     void add_2d_gaussian_beam_with_emittance( uint32_t N, double I, double q, double m,
00965                                               double a, double b, double e,
00966                                               double Ex, double x0 );
00967 };
00968 
00969 
00970 
00982 class ParticleDataBase3D : public ParticleDataBasePP<ParticleP3D> {
00983 
00984 
00985 public:
00986 
00987 /* ************************************** *
00988  * Constructors and destructor            *
00989  * ************************************** */
00990 
00993     ParticleDataBase3D() {}
00994 
00997     ~ParticleDataBase3D() {}
00998 
00999 /* ************************************** *
01000  * Particle beam definition               *
01001  * ************************************** */
01002 
01026     void add_cylindrical_beam_with_energy( uint32_t N, double J, double q, double m, 
01027                                            double E, double Tp, double Tt, Vec3D c, 
01028                                            Vec3D dir1, Vec3D dir2, double r );
01029 
01052     void add_cylindrical_beam_with_velocity( uint32_t N, double J, double q, double m, 
01053                                              double v, double dvp, double dvt, Vec3D c, 
01054                                              Vec3D dir1, Vec3D dir2, double r );
01055 
01064     void add_rectangular_beam_with_energy( uint32_t N, double J, double q, double m, 
01065                                            double E, double Tp, double Tt, Vec3D c, 
01066                                            Vec3D dir1, Vec3D dir2, double size1, double size2 );
01067 
01076     void add_rectangular_beam_with_velocity( uint32_t N, double J, double q, double m, 
01077                                              double v, double dvp, double dvt, Vec3D c, 
01078                                              Vec3D dir1, Vec3D dir2, double size1, double size2 );
01079 
01080 
01098     void add_3d_KV_beam_with_emittance( uint32_t N, double I, double q, double m,
01099                                         double ay, double by, double ey,
01100                                         double az, double bz, double ez,
01101                                         double Ex, double x0, double y0, double z0 );
01102 
01120     void add_3d_gaussian_beam_with_emittance( uint32_t N, double I, double q, double m,
01121                                               double ay, double by, double ey,
01122                                               double az, double bz, double ez,
01123                                               double Ex, double x0, double y0, double z0 );
01124 };
01125 
01126 
01127 
01128 
01129 
01130 
01131 #endif
01132 
01133 
01134 
01135 
01136 
01137 
01138 
01139 
01140 
01141 
01142 
01143 
01144 
01145 
01146 
01147 
01148 
01149 

Generated on Thu Apr 21 13:39:21 2011 for IBSimu by  doxygen 1.4.7