_kod źródłowy :

// prosta klasa implementująca ładunek
class el_charge {
public:
 double q;     -- ładunek
        vtx3d  pos;  - pozycja
        el_charge(const vtx3d&
_pos,double _q = -1.0)
         : q(_q)
                , pos(_pos) {}
        el_charge(float x,float y,float
z,double _q = -1.0)
         : q(_q)
                {
                   pos.set(x,y,z);
                }
};

void calc_field(const vtx3d& pos, vtx3d&
res);
float cos_between_vtx3d(const vtx3d&
a,const vtx3d&b);
void calc_path(el_charge& c,vtx3d& dir);

// tablica ładunków
el_charge* el_ch_t[20] = { 0 };

int el_ch_count = 0;


int dir_count = 14;
void get_dir(int dir,vtx3d& res);
// tu ciekawostka.
Z kazdego obiektu "wypuszczamy" linie
pola w <dir_count>
kierunkach, get_dir zwraca nam taki
wektor w i-tym kierunku,
niestety nie mam procedury dzialacej
kule na n czesci. Jak
masz jakis pomysl, to podeslij mi
procedurke get_dir ktora zwarca
n roznych kierunkow rozlozonych
rownomiernie po calej powierzchni
kuli.

void calc_path(el_charge& c,vtx3d& dir);

int fds = 0;

void calc_model()
{
  if( !fds ) {
             fds = 1

x     y     z   q
             el_ch_t[0] = new
el_charge(-0.5,-0.5 , 0,-5);
             el_ch_t[1] = new
l_charge( 0.0, 0.5 , 0,-3);
             el_ch_t[2] = new
el_charge(+0.5,-0.5 , 0,-1);

             el_ch_count = 3;

        }

        vtx3d dir;
        int i,j;
        el_charge* e;
        for( i = 0;i < el_ch_count;
i++ ) { // dla kazdego ladunku
             e = el_ch_t[i];
             for( j = 0; j < dir_count;
j++ ) {  // dla kazdego kierunku
                  get_dir(j,dir);    //
pobierz i-ty kierunek

                  calc_path(*e,dir);
// policz linie pola z *e ladunku w kierunku dir
             }
        }
}

// to jak pisalem: niezbyt logiczna na
twardo wpisana procedura zwracajaca
kierunki

void get_dir(int dir,vtx3d& res)
{
    switch( dir ) {
    case 0:
         res.set(-1.0, 0.0, 0.0);
         break;
    case 1:
         res.set( 1.0, 0.0, 0.0);
         break;
    case 2:
         res.set( 0.0,-1.0, 0.0);
         break;
    case 3:
         res.set( 0.0, 1.0, 0.0);
         break;
    case 4:
         res.set( 0.0, 0.0,-1.0);
         break;
    case 5:
         res.set( 0.0, 0.0, 1.0);
         break;

    case 6:
         res.set(-1.0,-1.0,-1.0);
         break;
    case 7:
         res.set( 1.0,-1.0,-1.0);
         break;
    case 8:
         res.set( 1.0, 1.0,-1.0);
         break;
    case 9:
         res.set(-1.0, 1.0,-1.0);
         break;
    case 10:
         res.set(-1.0,-1.0, 1.0);
         break;
    case 11:
         res.set( 1.0,-1.0, 1.0);
         break;
    case 12:
         res.set( 1.0, 1.0, 1.0);
         break;
    case 13:
         res.set(-1.0, 1.0, 1.0);
         break;
    }
}

//ilosc krokow rozniczkowania na jednostke
#define STEP_LEN  100.0
#define STEP (1.0/STEP_LEN)


// ze wzoru cosinusow
  c^2 = a^2 + b^2 - 2abcos(alpha)
                     (c^2 - a^2 - b^2)
 cos alpha =  ---------------------
                         -2 a b
float cos_betwe en_vtx3d(const vtx3d&
a,const vtx3d&b)
{
    vtx3d c   = a;
          c  -= b;
    float C   = c.length_sqr();    //
length_sqr() zwraca dlugosc wektora do
kwadratu
    float A   = a.length_sqr();
    float B   = b.length_sqr();

    return ( C - A - B ) /
            - 2 * (float)sqrt( A*B);
}

void calc_path(el_charge& c,vtx3d& dir)
{
    vtx3d pos = c.pos, cpos;
    vtx3d pdir, *va,*vb,v;
    float q = (float)c.q;
    bool end = false;
    dir *= (float)STEP;
    cpos = pos;
    pos += dir;
    int i = 0;
    while( !end ) {
// zapamietaj poprzedni kierunek

         if( i++ % 4 == 0 ) {
// tworz widzialna krawedz co 4 pomiary pola
// ta procedurka wpisuje do renderowanego "swiata" krawedz
edge_add(

vtx_add(pos.x,pos.y,pos.z),

vtx_add(cpos.x,cpos.y,cpos.z),
0
);
               // zapamietaj poprzednia pozycje
                cpos = pos;
         }

// policz wektor pola w punkcie pos do zmiennej dir
         calc_field(pos,dir);
         if( q > 0.0 )
             dir.invert();    //
v.invert() <==>  v * = -1;
         v = dir;              // zapamietaj wektor pola
         dir /= dir.length();    //wyciagnij wektor jednostkowy z dir
// o dlugosci jeden teraz |dir| = 1
         dir *= (float)STEP; // pomnoz przez nasze STEP rozniczkowania

             // jesli kat pomiedzy dir a pdir > 90 to znaczy, ze przeszlismy wlasnie przez ladunek i wektor pola odwrocil kierunek
         if( cos_between_vtx3d(dir,pdir)
< 0.0 ||
// lub jesli pozycja wyszla poza nasz ulubiony szescian
fabs(pos.x) >= 1.0 ||
fabs(pos.y) >= 1.0 ||
fabs(pos.z) >= 1.0 ) {
// konczymy liczenie sciezki
             end = true;
         } else /* q < 0.0 */   {
// w przeciwnym wypadku przesuwamy sie o wektor dir
pos += dir;
         }
    }
}

void calc_field(const vtx3d& pos, vtx3d&
res)
{
    vtx3d r,v;
    float rl;
    el_charge* ei;
    // ustawiamy rezulatat na zero
    res.set(0,0,0);
    // i superponujemy go ze potencjałami pola
       //pochodzacymi od kazdego ladunku w naszym ukladzie
    for( int i = 0; i < el_ch_count;
i++ ) if( ( ei = el_ch_t[i]) ){
        r  = pos;
        r -= ei->pos;
        rl = r.length();
        v  = r;
        assert(rl != 0.0);
        v /= rl;
        v *= 100.0*ei->q / (rl*rl);
        res -= v;
    }
}

copyright 2002