initialization.m

带电粒子在电磁场中运动轨迹模拟程序

function initialization( handles )
% initializes the vr world


%     Copyright (C) 2007 
%
%       mag. David Erzen
%       Faculty of Mechanical Engineering
%       LECAD Laboratory
%       Askerceva 6
%       1000 Ljubljana
%       SLOVENIA
%       contact email: [e-mail]
%       
%       Prof. John P. Verboncoeur
%       Plasma Theory and Simulation Group
%       University of California
%       Berkeley, CA 94720-1730 USA
%       
% 
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
% 
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%     GNU General Public License for more details.
% 
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <http://www.gnu.org/licenses/>.
    


    global_variables;

    % calculate E field and B field lines
    absolute_time = -1;
    [Bfield_lines Efield_lines] = BE_splines( handles );
    
    % set virtual world
    make_vr( 'sim1.wrl', Bfield_lines, Efield_lines );
    % TODO: the name of the virtual world file can be optional
    
    % set trajectory and guide center memory initial points
    trajectory = particle_initPos;
    clear global guide_center;    
    traj_set = get( handles.rbuttTraject, 'Value' ); 
    guide_set = get( handles.rbuttGuideCenter, 'Value' );
    vectorForce_set = get( handles.rbuttonVectorForce, 'Value' );
    vectorVel_set = get( handles.rbuttonVectorVelocity, 'Value' );
end

% ********** calculate magnetic and electric field lines ******************
function [Bfield_lines Efield_lines] = BE_splines( handles )

    global_variables;

    % initialize magnetic field lines
    B_parameters = str2num( get( handles.edtMagField,'String' ) );    
        
    switch get( handles.popMagField,'Value' )
        case BFIELD_NONE                                            % no magnetic field is selected
            spline_pos = [ 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2; 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2 ]; % starting point of splines, which represent field lines                
            Bline_length = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ];

        case BFIELD_CONSTANT                                        % constant magnetic filed in z direction
            spline_pos = [ 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2; 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2 ]; % starting point of splines, which represent field lines    
            Bline_length = [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ];

        case BFIELD_EXPONENT                                        % magnetic field in z direction, exponential fall in x direction
            spline_pos = [ 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2; 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2 ]; % starting point of splines, which represent field lines    
            Bline_length = [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ];
        
        case BFIELD_TORUS                                           % magnetic field in a torus           
            spline_pos = zeros(18, 3);
            temp_radius = zeros(18,1);
            for i=1:9
                spline_pos(i,1:3) = [B_parameters(2) + B_parameters(3)*cos( 2*pi*i/9 ) , 0,   B_parameters(3)*sin( 2*pi*i/9 ) ];
                temp_radius(i,1) = sqrt(spline_pos(i,1)^2 + spline_pos(i,3)^2);
            end            
            spline_pos(10:18, 1:3) = spline_pos(1:9, 1:3);  
            spline_pos(10:18, 1) = -spline_pos(1:9, 1);     % made for
            %calculating magnetic field lines in opposite direction
            temp_radius(10:18, 1) = temp_radius(1:9,1);
            
            Bline_length = 40*temp_radius;
        
        case BFIELD_TORUS_POLOIDAL                                  % magnetic field in a torus with poloidal field
            for i=1:9
                spline_pos(i,1:3) = [B_parameters(2) - B_parameters(4)*cos( pi/2+pi*i/9 ) , 0,   B_parameters(4)*sin( pi/2+pi*i/9 ) ];
                % temporary positions:
                spline_pos(i+9,1:3) = [-B_parameters(2) - B_parameters(4)*cos( -pi/2 + pi*i/9 ) , 0,   B_parameters(4)*sin( -pi/2 + pi*i/9 ) ];
            end
            % spline_pos(10:18, 1:3) = spline_pos(1:9, 1:3); 
            % spline_pos(10:18, 1) = -spline_pos(1:9, 1);     % made for
            % calculating magnetic field lines in opposite direction
            Bline_length = 1.6*(B_parameters(2) +B_parameters(4) )* [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ];
           
        case BFIELD_CONCENTRIC
            spline_pos = [ B_parameters(2) 0 0.1; B_parameters(2) 0 1; B_parameters(2) 0 2; B_parameters(2) 0 3; B_parameters(2) 0 4; B_parameters(2) 0 -1; B_parameters(2) 0 -2; B_parameters(2) 0 -3; B_parameters(2) 0 -4; 3*B_parameters(2) 0 0.1; 3*B_parameters(2) 0 1; 3*B_parameters(2) 0 2; 3*B_parameters(2) 0 3; 3*B_parameters(2) 0 4; 3*B_parameters(2) 0 -1; 3*B_parameters(2) 0 -2; 3*B_parameters(2) 0 -3; 3*B_parameters(2) 0 -4]; % starting point of splines, which represent field lines    
            Bline_length = 12*B_parameters(2)*[100, 100, 100, 100, 100, 100, 100, 100, 100, 300, 300, 300, 300, 300, 300, 300, 300, 300 ];
            
        case BFIELD_CUSTOM                                          % a user defined magnetic field is selected 
            %spline_pos = [ 0.01 0.01 3.01; 0.01 -0.01 3.01; -0.01 -0.01 3.01; -0.01 0.01 3.01; 0.01 0.01 2.99; 0.01 -0.01 2.99; -0.01 0.01 2.99; -0.01 -0.01 2.99; 0.01 -0.01 2.99; 0.01 0.01 3; -0.01 0.01 3; -0.01 -0.01 3; 0 -0.01 2.99; 0 0.01 2.99 ; 0 0.01 3.01; -0.01 -0.01 3; 0 -0.01 2.99; 0 0 3.01 ]; % starting point of splines, which represent field lines            
            %spline_pos = [ 0.001 0.001 3.01; 0.001 -0.001 3.01; -0.001 -0.001 3.01; -0.001 0.001 3.01; 0.001 0.001 7.99; 0.001 -0.001 7.99; -0.001 0.001 7.99; -0.001 -0.001 7.99; 0.001 -0.001 7.99; 0.001 0.001 3.01; -0.001 0.001 3.01; -0.001 -0.001 3.01; 0 -0.001 7.99; 0 0.001 7.99 ; 0 0.001 3.01; -0.001 -0.001 3; 0 -0.001 7.99; 0.001 0 3.01 ]; % starting point of splines, which represent field lines            
            %spline_pos = [ 1.0015 1.0015 0.01; 1.0015 0.9985 0.01; 0.9985 0.9985 0.01; 0.9985 1.0015 0.01; 1.0015 1.0015 4.99; 1.0015 0.9985 4.99; 0.9985 1.0015 4.99; 0.9985 0.9985 4.99; 1.0015 0.9985 4.99; 1.0015 1.0015 0.01; 0.9985 1.0015 0.01; 0.9985 0.9985 0.01; 1 0.9985 4.99; 1 1.0015 4.99 ; 1 1.0015 0.01; 1.0 1.0 8; 1 0.9985 4.99; 1.0015 1 0.01 ];
            %spline_pos = [ 1 1.001 0.008; 1.0 0.999 0.008; 1.0 1.0 0.008; 1.001 1.001 0.008; 1.001 0.999 0.008; 0.999 0.999 0.008; 0.999 1.001 0.008; 1.001 1.0 0.008; 0.999 1.0 0.008; 1.005 1.005 0.008; 1.005 0.995 0.008; 0.995 1.005 0.008; 0.995 0.995 0.008; 1.005 1.0 0.008; 0.995 1.0 0.008; 1 0.995 0.008; 1 1.005 0.008 ; 1.0 1.0 5.995]; % starting point of splines, which represent field lines            
            spline_pos = [ 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2; 0 0 0.1; 1 0 -1; 1 1 0; 0 1 1; -1 1 2; -1 0 1; -1 -1 0; 0 -1 -1; 1 -1 -2 ]; % starting point of splines, which represent field lines  
            %spline_pos = [ 0.1 0.1 3.1; 0.1 -0.1 3.1; -0.1 -0.1 3.1; -0.1 0.1 3.1; 0.1 0.1 2.9; 0.1 -0.1 2.9; -0.1 0.1 2.9; -0.1 -0.1 2.9; 0.1 -0.1 2.9; 0.1 0.1 3; -0.1 0.1 3; -0.1 -0.1 3; 0 -0.1 2.9; 0 0.1 2.9 ; 0 0.1 3.1; -0.1 -0.1 3; 0 -0.1 2.9; 0 0 3.1 ]; % starting point of splines, which represent field lines            
            Bline_length = [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ];
            
        case BFIELD_MIRROR
            dl = 0.1;
            dz = 0.3;
            spline_pos = [ 0, 0, B_parameters(1)/2-dz; dl, 0, B_parameters(1)/2-dz; 0, dl, B_parameters(1)/2-dz; -dl, 0, B_parameters(1)/2-dz; 0, -dl, B_parameters(1)/2-dz; 0, 0, -B_parameters(1)/2+dz; dl, 0, -B_parameters(1)/2+dz; 0, dl, -B_parameters(1)/2+dz; -dl, 0, -B_parameters(1)/2+dz; 0, -dl, -B_parameters(1)/2+dz; dl, dl, B_parameters(1)/2-dz; -dl, dl, B_parameters(1)/2-dz; -dl, -dl, B_parameters(1)/2-dz; dl, -dl, B_parameters(1)/2-dz; dl, dl, -B_parameters(1)/2+dz; -dl, dl, -B_parameters(1)/2+dz; -dl, -dl, -B_parameters(1)/2+dz; dl, -dl, -B_parameters(1)/2+dz ];
            
            if B_parameters(2) > 0  % only if moonpole is positive a field line should begin near by otherwise it will collapse into the monopole
                top_length = 1;
            else
                top_length = 1e-6;
            end;
            
            if B_parameters(3) > 0 
                bottom_length = 1;
            else
                bottom_length = 1e-6;
            end;