Потеря памяти в программе

@Default · Регистрация: 08.04.2012

Author24 — интернет-сервис помощи студентам

Подскажите пожалуйста, как отловить в каком месте происходит потеря памяти? Программа считает, но если открыть диспетчер задач, то файл подкачки растет до 2ГБ и программа перестает считать, а там внутри цикл, который должен работать долго.
-------------------------------
Приложен txt с кодом. Коду не пугаться =)

Croessmah · 24.10.2012, 18:41

Сообщение от Default

Как я догадываюсь основная проблема у меня в том, что когда я из функции get_B() вызываю функцию direct_multiple() у меня происходит накопление памяти. Как от него избавиться не ясно. Когда перед return прописываю delete[] matrfi и matrthfi, в которые я передаю direct_multiple() вылазит фатал, из-за чего конкретно и в чем косяк разобраться пока что не выходит.

Исходя из сказанного Вами, могу сделать вывод - Вы просто не понимаете как работать с памятью. Посему еще раз говорю - откройте книжку и почитайте. Можете взять любую книгу по C++ и почитать главу "Динамическая память(куча)"

@Default · 24.10.2012, 18:43 **[ТС]**

Не понимаю, сейчас буду смотреть видео лекцию / искать книги.

alsav22 · 25.10.2012, 04:53

Добавлено через 4 часа 28 минут
Решение:

Код

C++

#include <math.h>
#include <time.h>
#include <stdio.h>
 
double pi = 3.14159265358979323846264338327;
long N1rand=101;
long N2rand=302;
 
//double* matrfi = new double[4];
//double* matrthfi = new double[4];
 
double* new_matrixS = new double[4];
double* matr = new double[4];   
double* Ef = new double[4];
double* g = new double[2];
 
struct work
{
    double answer[2];
    double dispersion[2];
    double srednee;
    int max_dlina;
};
 
double rand(void);
double* direct_multiple(double angle1, double angle2);
double* funkg(double koord[],double mu, double lambda);
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[]);
//изменил название на matr2, чтобы не совпадало с указателем
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda);
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda);
 
int main()
{
    int ispitania1  = 1000000;
    int ispitania10 = 10000000;
    work s4et;
    double* x0 = new double[2];
    double* x = new double[2];
    double* x1 = new double[2];
    double* z = new double[2];
    double R0 = 1;
    double R1 = 1;
    clock_t time;
    double lambda;
    double mu;
    double alpha;
 
    x0[0]=0;
    x0[1]=0;
 
    x1[0]=1;
    x1[1]=1;
 
    mu=1;
    alpha=1;
    lambda=mu*(alpha-1);
    for(int i=0; i<25;i++)
    {
        x[0]= -0.7+(0.1*i);
        x[1]= -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        //double* temp;
        z=funkg(x,mu,lambda);
        //temp = z;
        time = clock();
        s4et = tree(x,x0,x1,R0,R1,ispitania1,mu,lambda);
        time = clock() - time;
        //fprintf(lame2_1_1,"%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",i+1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
        //if (abs(s4et.answer[0]-z[0])<s4et.dispersion[0]){printf("it is ok");}
        //else {printf("it is not ok");}
        printf("%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
    }
    
    delete [] x0;
    delete [] x;
    delete [] x1;
    if (z == g) delete [] z;
    else 
        {
            delete [] z;
            delete [] g;
        }
    
    if (matr == Ef) 
    {
        delete [] matr; 
        delete [] new_matrixS;
    }
    else 
        if (matr == new_matrixS) 
        {
            delete [] matr;
            delete [] Ef;
        }
        else
        {
            delete [] matr;
            delete [] Ef;
            delete [] new_matrixS;
        }
   
 
    getchar();
    return 0;
}
 
double rand(void)/* the initial data(type long):  N1rand=101; N2rand=302; */
{ 
  double x;
  long k;
  k=x=((double)( N1rand=(( N1rand*=1373)%1000919))/1000919.0+
  (double)( N2rand=(( N2rand*=1528)%1400159))/1400159.0 );
  return(x-k);
}
 
double* direct_multiple(double angle1, double angle2)
{
    //double* new_matrixS = new double[4];// пока что с двойкой, потом двойка будет входить в переменную для весов
    
    new_matrixS[0]=2*cos(angle1)*cos(angle2);
    new_matrixS[1]=2*cos(angle1)*sin(angle2);
    new_matrixS[2]=2*sin(angle1)*cos(angle2);
    new_matrixS[3]=2*sin(angle1)*sin(angle2);
    return new_matrixS;
}
 
 
 
double* funkg(double koord[],double mu, double lambda)//,double Z, double fi)
{
    
    //double* g = new double[2];
    g[0] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[0],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    g[1] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[1],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    return g;
}
 
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[])
{
    
    double thetta;
    thetta = atan2(y[1]-Rc[1],y[0]- Rc[0]);
 
    double* matrfi = NULL;// = new double[4];
    double* matrthfi = NULL;// = new double[4];
    
    //double* matr = new double[4]; 
    //double* Ef = new double[4];
    
    Ef[0]=1;
    Ef[1]=0;
    Ef[2]=0;
    Ef[3]=1;
    
    double p0 = 0;
 
    p0 = (alpha*Z - R)/R;
 
    double a1; // случайная величина от 0 до 1 используется для вычисления матрицы B
    double a2; // случайная величина от 0 до 1 используется для вычисления матрицы B если a1 > alpha/(2+alpha)
 
    a1 = rand();
    a2 = rand();
 
    if (thetta < 0) {thetta = thetta + 2*pi;}
    
    //matrfi = direct_multiple(fi,fi);
    //matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
    
    if (a1 < alpha/(2+alpha))
    {
        matrfi = direct_multiple(fi,fi);
        matr = matrfi;
    }
    else if (2*R - alpha*Z > 0) 
         {
            if (a2 < (2*R - alpha*Z)/2*R)
            {
                matr = Ef;
            }
            else 
            {
                matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
                matr = matrthfi;
            }
         }
         else if (a2 < (alpha*Z - 2*R)/(2*R*p0))
              {
                    matr[0] = -1*p0;
                    matr[1] = 0;
                    matr[2] = 0;
                    matr[3] = -1*p0;
              }
              else
              {
                    matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
                    matr[0] = matrthfi[0]*p0;
                    matr[1] = matrthfi[1]*p0;
                    matr[2] = matrthfi[2]*p0;
                    matr[3] = matrthfi[3]*p0;
              }
    
    return matr;
}
 
 
 
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda)
{//t0 - na4alo koord, tn - to4ka reshenia, tr - slu4ainaia to4ka
    
    double* matr = new double[4];
        
    double costt = (tn[0]-t0[0]+Z*cos(fi))/radius;
 
    double sintt = (tn[1]-t0[1]+Z*sin(fi))/radius;
 
    double cosfi = cos(fi);
 
    double sinfi = sin(fi);
 
    matr[0]=(((1 + 2*mu/(lambda + mu))-1) + 2*cosfi*cosfi + (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr[1]=(2*cosfi*sinfi + (Z/radius)*(sintt*cosfi + costt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr[2]=matr[1];
    
    matr[3]=(((1 + 2*mu/(lambda + mu))-1) + 2*sinfi*sinfi - (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    return matr;
}
 
 
 
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda)
{
    work resh;
    int k;
    int k2;
    double* t4 = new double[2]; //peremennaia to4ka dlia blujdania
    double* y = new double[2];
    int* obriv = new int[100];
    int obr;
    int max_obriv = 0;
    double* d = NULL;// = new double[2];
    double* matrix = NULL;//new double[4];
    double* E = new double[4];
    double* Q = new double[4];
    double* Qt = new double[4];
    double p = 0;
    double modz = 0;
    double N2 = N;
    double fp = 0;
    double srednee = 0;
    int oi = 0;
 
 
    double p0 = 0; // надо в том случае, если 2R-alpha*Z < 0 
 
    double fi;
    double fipi;
    double J;
    double psi;
    double Z;
    double Zs;
    double a; // alpha случайная величина от 0 до 1
    
    double fis; // угол фи штрих
    double thetta; // угол тетта
 
    double alpha = (mu + lambda)/mu;
 
        
    resh.answer[0]=0;
    resh.answer[1]=0;
    resh.dispersion[0]=0;
    resh.dispersion[1]=0;
 
    E[0]=1;
    E[1]=0;
    E[2]=0;
    E[3]=1;
 
    for (int pp=0; pp<100; pp++)
        obriv[pp]=0;
 
    if (sqrt(pow((tn[0]-t0[0]),2) + pow((tn[1]-t0[1]),2)) < R0) k2=0;
    else k2=1;
    for(int i=0; i<N; i++)
    {
        obr=0;
        k=k2;
        for(int s=0;s<4;s++)
        {
            Q[s]=E[s];
        }
        for(int s=0;s<2;s++)
        {
            t4[s]=tn[s];
        }
        while(k!=2)
        {
            if (k==0)
            {   
                fis = atan2(t4[1]-t0[1],t4[0]- t0[0]);
                if (fis < 0) {fis = fis + 2*pi;}
                p = sqrt(pow(t4[0]-t0[0],2)+ pow(t4[1]-t0[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R0,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R0,t0,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
                    
                    matrix = get_B(alpha,fipi,R0,t0,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t1[0]),2) + pow((y[1]-t1[1]),2)) < R1) 
                {
                    k=1;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    //Q=Qt*matr;
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //delete[] d;
                    //break;
                }
                
            } // if (k == 0)
            if (k==1)
            {
                fis = atan2(t4[1]-t1[1],t4[0]- t1[0]);
                p = sqrt(pow(t4[0]-t1[0],2)+ pow(t4[1]-t1[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R1,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R1,t1,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
 
                    matrix = get_B(alpha,fipi,R1,t1,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t0[0]),2) + pow((y[1]-t0[1]),2)) < R0)
                {
                    k=0;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //delete[] d;
                    
                    //break;
                }
            } //if (k ==1)
                
        } // while
    if (obr < 100) obriv[obr]++;
    }
    
    int obr2 = 99;
    while(obriv[obr2]==0) obr2--;
    resh.max_dlina = obr2+1;    
 
    for (int i=0; i<100; i++)
    {
        srednee+=(i+1)*obriv[i];
    }
    resh.srednee = srednee/N2;
        
    resh.answer[0]=resh.answer[0]/N;
    resh.answer[1]=resh.answer[1]/N;
    
    resh.dispersion[0]= 3*sqrt(resh.dispersion[0]/N - pow(resh.answer[0],2)) / sqrt(N2);
    resh.dispersion[1]= 3*sqrt(resh.dispersion[1]/N - pow(resh.answer[1],2)) / sqrt(N2);
 
    delete [] E;
    delete [] Q;
    delete [] Qt;
 
    delete [] t4;
    delete [] y;
    delete [] obriv;
    
    return resh;
}

@Default · 25.10.2012, 09:25 **[ТС]**

Огромное спасибо, теперь буду вникать и пытаться разбираться как оно теперь работает =)

alsav22 · 25.10.2012, 19:40

Сообщение от Default

Огромное спасибо, теперь буду вникать и пытаться разбираться как оно теперь работает =)

Работает оно, в общем, так же. Основное, что я сделал, это вынес выделение памяти из функций. Особенно это нужно было для new_matrixS и matr (выделение памяти под них, в основном, и пожирало память). Что ещё при этом пришлось изменить, так это, в get_B() присвоение matr указателей на память, заполненную в direct_multiple(). Там шло подряд два вызова direct_multiple(). Естественно, первое содержимое памяти, после второго вызова, уже нельзя было использовать. Но в коде, применение этой памяти разнесено по if else, поэтому и сами вызовы можно разнести. Вот тут код, в котором я сделал небольшие изменения. В коде, функциия matr2() (бывшая matr()), насколько я понял, не используется, но если будет, то выделение памяти из неё тоже нужно вынести (что я и сделал). Когда занимался этим кодом, то пришла мысль, что, во многих местах, применение динамической памяти совсем не обязательно. Вполне можно обойтись и статической.

Код

C++

#include <math.h>
#include <time.h>
#include <stdio.h>
 
double pi = 3.14159265358979323846264338327;
long N1rand=101;
long N2rand=302;
 
double* new_matrixS = new double[4]; // указатель, который используется в direct_multiple()
double* matr = new double[4]; // указатель, который используется в get_B()
double* del_matr = matr; // указатель для освобождения памяти под matr     
double* Ef = new double[4]; // указатель, который используется в get_B()
double* g = new double[2]; // указатель, который используется в funkg()
 
double* matr3 = new double[4]; // указатель, который используется в matr2()
 
struct work
{
    double answer[2];
    double dispersion[2];
    double srednee;
    int max_dlina;
};
 
double rand(void);
double* direct_multiple(double angle1, double angle2);
double* funkg(double koord[],double mu, double lambda);
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[]);
//изменил имя на matr2, чтобы не совпадало с указателем
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda);
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda);
 
int main()
{
    int ispitania1  = 1000000;
    int ispitania10 = 10000000;
    work s4et;
    double* x0 = new double[2];
    double* x = new double[2];
    double* x1 = new double[2];
    double* z = NULL; /*= new double[2]*/;
    double R0 = 1;
    double R1 = 1;
    clock_t time;
    double lambda;
    double mu;
    double alpha;
 
    x0[0]=0;
    x0[1]=0;
 
    x1[0]=1;
    x1[1]=1;
 
    mu=1;
    alpha=1;
    lambda=mu*(alpha-1);
    for(int i=0; i<25;i++)
    {
        x[0]= -0.7+(0.1*i);
        x[1]= -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        //double* temp;
        z=funkg(x,mu,lambda);
        //temp = z;
        time = clock();
        s4et = tree(x,x0,x1,R0,R1,ispitania1,mu,lambda);
        time = clock() - time;
        //fprintf(lame2_1_1,"%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",i+1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
        //if (abs(s4et.answer[0]-z[0])<s4et.dispersion[0]){printf("it is ok");}
        //else {printf("it is not ok");}
        printf("%d %4.1f %8.5f %8.5f %8.5f %8.5f %8.5f %8.5f %d %4.1f\n",1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
    }
    
    delete [] x0;
    delete [] x;
    delete [] x1;
    
    delete [] del_matr;
    delete [] Ef;   
    delete [] new_matrixS;
    delete [] g;
    delete [] matr3;
            
 
    getchar();
    return 0;
}
 
double rand(void)/* the initial data(type long):  N1rand=101; N2rand=302; */
{ 
  double x;
  long k;
  k=x=((double)( N1rand=(( N1rand*=1373)%1000919))/1000919.0+
  (double)( N2rand=(( N2rand*=1528)%1400159))/1400159.0 );
  return(x-k);
}
 
double* direct_multiple(double angle1, double angle2)
{
    // пока что с двойкой, потом двойка будет входить в переменную для весов
    
    new_matrixS[0]=2*cos(angle1)*cos(angle2);
    new_matrixS[1]=2*cos(angle1)*sin(angle2);
    new_matrixS[2]=2*sin(angle1)*cos(angle2);
    new_matrixS[3]=2*sin(angle1)*sin(angle2);
    return new_matrixS;
}
 
 
double* funkg(double koord[],double mu, double lambda)//,double Z, double fi)
{
    g[0] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[0],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    g[1] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[1],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    return g;
}
 
 
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[])
{
    double thetta;
    thetta = atan2(y[1]-Rc[1],y[0]- Rc[0]);
 
    double p0 = 0;
 
    p0 = (alpha*Z - R)/R;
 
    double a1; // случайная величина от 0 до 1 используется для вычисления матрицы B
    double a2; // случайная величина от 0 до 1 используется для вычисления матрицы B если a1 > alpha/(2+alpha)
 
    a1 = rand();
    a2 = rand();
 
    if (thetta < 0) thetta = thetta + 2*pi;
    
    if (a1 < alpha/(2+alpha)) matr = direct_multiple(fi,fi);
    else if (2*R - alpha*Z > 0) 
         {
            if (a2 < (2*R - alpha*Z)/2*R)
            {
                Ef[0]=1;
                Ef[1]=0;
                Ef[2]=0;
                Ef[3]=1;
                
                matr = Ef;
            }
            else matr = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
         }
         else if (a2 < (alpha*Z - 2*R)/(2*R*p0))
              {
                    matr[0] = -1*p0;
                    matr[1] = 0;
                    matr[2] = 0;
                    matr[3] = -1*p0;
              }
              else
              {
                    double* matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
                    
                    matr[0] = matrthfi[0]*p0;
                    matr[1] = matrthfi[1]*p0;
                    matr[2] = matrthfi[2]*p0;
                    matr[3] = matrthfi[3]*p0;
              }
    
    return matr;
}
 
 
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda)
{//t0 - na4alo koord, tn - to4ka reshenia, tr - slu4ainaia to4ka
    
    //double* matr3 = new double[4];
        
    double costt = (tn[0]-t0[0]+Z*cos(fi))/radius;
 
    double sintt = (tn[1]-t0[1]+Z*sin(fi))/radius;
 
    double cosfi = cos(fi);
 
    double sinfi = sin(fi);
 
    matr3[0]=(((1 + 2*mu/(lambda + mu))-1) + 2*cosfi*cosfi + (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[1]=(2*cosfi*sinfi + (Z/radius)*(sintt*cosfi + costt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[2]=matr[1];
    
    matr3[3]=(((1 + 2*mu/(lambda + mu))-1) + 2*sinfi*sinfi - (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    return matr3;
}
 
 
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda)
{
    work resh;
    int k;
    int k2;
    double* t4 = new double[2]; //peremennaia to4ka dlia blujdania
    double* y = new double[2];
    int* obriv = new int[100];
    int obr;
    int max_obriv = 0;
    double* d = NULL;
    double* matrix = NULL;
    double* E = new double[4];
    double* Q = new double[4];
    double* Qt = new double[4];
    double p = 0;
    double modz = 0;
    double N2 = N;
    double fp = 0;
    double srednee = 0;
    int oi = 0;
 
    double p0 = 0; // надо в том случае, если 2R-alpha*Z < 0 
 
    double fi;
    double fipi;
    double J;
    double psi;
    double Z;
    double Zs;
    double a; // alpha случайная величина от 0 до 1
    
    double fis; // угол фи штрих
    double thetta; // угол тетта
 
    double alpha = (mu + lambda)/mu;
 
    resh.answer[0]=0;
    resh.answer[1]=0;
    resh.dispersion[0]=0;
    resh.dispersion[1]=0;
 
    E[0]=1;
    E[1]=0;
    E[2]=0;
    E[3]=1;
 
    for (int pp=0; pp<100; pp++)
        obriv[pp]=0;
 
    if (sqrt(pow((tn[0]-t0[0]),2) + pow((tn[1]-t0[1]),2)) < R0) k2=0;
    else k2=1;
    for(int i=0; i<N; i++)
    {
        obr=0;
        k=k2;
        for(int s=0;s<4;s++)
        {
            Q[s]=E[s];
        }
        for(int s=0;s<2;s++)
        {
            t4[s]=tn[s];
        }
        while(k!=2)
        {
            if (k==0)
            {   
                fis = atan2(t4[1]-t0[1],t4[0]- t0[0]);
                if (fis < 0) {fis = fis + 2*pi;}
                p = sqrt(pow(t4[0]-t0[0],2)+ pow(t4[1]-t0[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R0,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R0,t0,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
                    
                    matrix = get_B(alpha,fipi,R0,t0,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t1[0]),2) + pow((y[1]-t1[1]),2)) < R1) 
                {
                    k=1;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    //Q=Qt*matr;
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //delete[] d;
                    //break;
                }
                
            } // if (k == 0)
            if (k==1)
            {
                fis = atan2(t4[1]-t1[1],t4[0]- t1[0]);
                p = sqrt(pow(t4[0]-t1[0],2)+ pow(t4[1]-t1[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R1,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R1,t1,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
 
                    matrix = get_B(alpha,fipi,R1,t1,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t0[0]),2) + pow((y[1]-t0[1]),2)) < R0)
                {
                    k=0;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //delete[] d;
                    
                    //break;
                }
            } //if (k ==1)
                
        } // while
    if (obr < 100) obriv[obr]++;
    }
    
    int obr2 = 99;
    while(obriv[obr2]==0) obr2--;
    resh.max_dlina = obr2+1;    
 
    for (int i=0; i<100; i++)
    {
        srednee+=(i+1)*obriv[i];
    }
    resh.srednee = srednee/N2;
        
    resh.answer[0]=resh.answer[0]/N;
    resh.answer[1]=resh.answer[1]/N;
    
    resh.dispersion[0]= 3*sqrt(resh.dispersion[0]/N - pow(resh.answer[0],2)) / sqrt(N2);
    resh.dispersion[1]= 3*sqrt(resh.dispersion[1]/N - pow(resh.answer[1],2)) / sqrt(N2);
 
    delete [] E;
    delete [] Q;
    delete [] Qt;
 
    delete [] t4;
    delete [] y;
    delete [] obriv;
    
    return resh;
}

Добавлено через 15 минут
Вот вариант, когда динамически память выделяется только под matr:

Код

C++

#include <math.h>
#include <time.h>
#include <stdio.h>
 
double pi = 3.14159265358979323846264338327;
long N1rand=101;
long N2rand=302;
 
 
double* matr = new double[4]; // указатель, который используется в get_B()
double* del_matr = matr; // указатель для освобождения памяти под matr 
 
double new_matrixS[4]; 
double Ef[4]; 
double g[2]; 
double matr3[4]; 
 
struct work
{
    double answer[2];
    double dispersion[2];
    double srednee;
    int max_dlina;
};
 
double rand(void);
double* direct_multiple(double angle1, double angle2);
double* funkg(double koord[],double mu, double lambda);
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[]);
//изменил имя на matr2, чтобы не совпадало с указателем
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda);
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda);
 
int main()
{
    int ispitania1  = 1000000;
    int ispitania10 = 10000000;
    work s4et;
    double x0[2];
    double x[2];
    double x1[2];
    double* z = NULL;
    double R0 = 1;
    double R1 = 1;
    clock_t time;
    double lambda;
    double mu;
    double alpha;
 
    x0[0]=0;
    x0[1]=0;
 
    x1[0]=1;
    x1[1]=1;
 
    mu=1;
    alpha=1;
    lambda=mu*(alpha-1);
    for(int i=0; i<25;i++)
    {
        x[0] = -0.7+(0.1*i);
        x[1] = -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        //double* temp;
        z = funkg(x, mu, lambda);
        //temp = z;
        time = clock();
        s4et = tree(x, x0, x1, R0, R1, ispitania1, mu, lambda);
        time = clock() - time;
        //fprintf(lame2_1_1,"%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",i+1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
        //if (abs(s4et.answer[0]-z[0])<s4et.dispersion[0]){printf("it is ok");}
        //else {printf("it is not ok");}
        printf("%d %4.1f %8.5f %8.5f %8.5f %8.5f %8.5f %8.5f %d %4.1f\n",1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
    }
    
    
    delete [] del_matr;
        
    getchar();
    return 0;
}
 
 
double rand(void)/* the initial data(type long):  N1rand=101; N2rand=302; */
{ 
  double x;
  long k;
  k=x=((double)( N1rand=(( N1rand*=1373)%1000919))/1000919.0+
  (double)( N2rand=(( N2rand*=1528)%1400159))/1400159.0 );
  return(x-k);
}
 
 
double* direct_multiple(double angle1, double angle2)
{
    // пока что с двойкой, потом двойка будет входить в переменную для весов
    
    new_matrixS[0]=2*cos(angle1)*cos(angle2);
    new_matrixS[1]=2*cos(angle1)*sin(angle2);
    new_matrixS[2]=2*sin(angle1)*cos(angle2);
    new_matrixS[3]=2*sin(angle1)*sin(angle2);
    return new_matrixS;
}
 
 
double* funkg(double koord[],double mu, double lambda)//,double Z, double fi)
{
    g[0] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[0],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    g[1] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[1],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    return g;
}
 
 
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[])
{
    double thetta;
    thetta = atan2(y[1]-Rc[1],y[0]- Rc[0]);
 
    double p0 = 0;
 
    p0 = (alpha*Z - R)/R;
 
    double a1; // случайная величина от 0 до 1 используется для вычисления матрицы B
    double a2; // случайная величина от 0 до 1 используется для вычисления матрицы B если a1 > alpha/(2+alpha)
 
    a1 = rand();
    a2 = rand();
 
    if (thetta < 0) thetta = thetta + 2*pi;
    
    if (a1 < alpha/(2+alpha)) matr = direct_multiple(fi,fi);
    else if (2*R - alpha*Z > 0) 
         {
            if (a2 < (2*R - alpha*Z)/2*R)
            {
                Ef[0]=1;
                Ef[1]=0;
                Ef[2]=0;
                Ef[3]=1;
                
                matr = Ef;
            }
            else matr = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
         }
         else if (a2 < (alpha*Z - 2*R)/(2*R*p0))
              {
                    matr[0] = -1*p0;
                    matr[1] = 0;
                    matr[2] = 0;
                    matr[3] = -1*p0;
              }
              else
              {
                    double* matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
                    
                    matr[0] = matrthfi[0]*p0;
                    matr[1] = matrthfi[1]*p0;
                    matr[2] = matrthfi[2]*p0;
                    matr[3] = matrthfi[3]*p0;
              }
    
    return matr;
}
 
 
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda)
{//t0 - na4alo koord, tn - to4ka reshenia, tr - slu4ainaia to4ka
    
    double costt = (tn[0]-t0[0]+Z*cos(fi))/radius;
 
    double sintt = (tn[1]-t0[1]+Z*sin(fi))/radius;
 
    double cosfi = cos(fi);
 
    double sinfi = sin(fi);
 
    matr3[0]=(((1 + 2*mu/(lambda + mu))-1) + 2*cosfi*cosfi + (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[1]=(2*cosfi*sinfi + (Z/radius)*(sintt*cosfi + costt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[2]=matr[1];
    
    matr3[3]=(((1 + 2*mu/(lambda + mu))-1) + 2*sinfi*sinfi - (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    return matr3;
}
 
 
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda)
{
    work resh;
    int k;
    int k2;
    int obr;
    int max_obriv = 0;
    
    double t4[2]; //peremennaia to4ka dlia blujdania
    double y[2];
    int obriv[100];
    
    double E[4];
    double Q[4];
    double Qt[4];
 
    double* d = NULL;
    double* matrix = NULL;
 
    double p = 0;
    double modz = 0;
    double N2 = N;
    double fp = 0;
    double srednee = 0;
    int oi = 0;
 
    double p0 = 0; // надо в том случае, если 2R-alpha*Z < 0 
 
    double fi;
    double fipi;
    double J;
    double psi;
    double Z;
    double Zs;
    double a; // alpha случайная величина от 0 до 1
    
    double fis; // угол фи штрих
    double thetta; // угол тетта
 
    double alpha = (mu + lambda)/mu;
 
    resh.answer[0]=0;
    resh.answer[1]=0;
    resh.dispersion[0]=0;
    resh.dispersion[1]=0;
 
    E[0]=1;
    E[1]=0;
    E[2]=0;
    E[3]=1;
 
    for (int pp=0; pp<100; pp++)
        obriv[pp]=0;
 
    if (sqrt(pow((tn[0]-t0[0]),2) + pow((tn[1]-t0[1]),2)) < R0) k2=0;
    else k2=1;
    for(int i=0; i<N; i++)
    {
        obr=0;
        k=k2;
        for(int s=0;s<4;s++)
        {
            Q[s]=E[s];
        }
        for(int s=0;s<2;s++)
        {
            t4[s]=tn[s];
        }
        while(k!=2)
        {
            if (k==0)
            {   
                fis = atan2(t4[1]-t0[1],t4[0]- t0[0]);
                if (fis < 0) {fis = fis + 2*pi;}
                p = sqrt(pow(t4[0]-t0[0],2)+ pow(t4[1]-t0[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R0,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R0,t0,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
                    
                    matrix = get_B(alpha,fipi,R0,t0,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t1[0]),2) + pow((y[1]-t1[1]),2)) < R1) 
                {
                    k=1;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    //Q=Qt*matr;
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    
                    //break;
                }
                
            } // if (k == 0)
            if (k==1)
            {
                fis = atan2(t4[1]-t1[1],t4[0]- t1[0]);
                p = sqrt(pow(t4[0]-t1[0],2)+ pow(t4[1]-t1[1],2));
 
                fi = rand()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R1,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R1,t1,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
 
                    matrix = get_B(alpha,fipi,R1,t1,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t0[0]),2) + pow((y[1]-t0[1]),2)) < R0)
                {
                    k=0;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    
                    //continue;
                }
                else
                {
                    d = funkg(y,mu,lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //break;
                }
            } //if (k ==1)
                
        } // while
    if (obr < 100) obriv[obr]++;
    }
    
    int obr2 = 99;
    while(obriv[obr2]==0) obr2--;
    resh.max_dlina = obr2+1;    
 
    for (int i=0; i<100; i++)
    {
        srednee+=(i+1)*obriv[i];
    }
    resh.srednee = srednee/N2;
        
    resh.answer[0]=resh.answer[0]/N;
    resh.answer[1]=resh.answer[1]/N;
    
    resh.dispersion[0]= 3*sqrt(resh.dispersion[0]/N - pow(resh.answer[0],2)) / sqrt(N2);
    resh.dispersion[1]= 3*sqrt(resh.dispersion[1]/N - pow(resh.answer[1],2)) / sqrt(N2);
 
    return resh;
}

@Default · 25.10.2012, 19:45 **[ТС]**

Появилась проблема, которую я сейчас пытаюсь отловить, в мэйне

z=funkg(x,mu,lambda);

time = clock();

s4et = tree(x,x0,x1,R0,R1,ispitania1,mu,lambda);

time = clock() - time;

после вызова функции tree, меняется значение z, что никак не допустимо для вычислений. Видимо z идет как указатель на funkg() и когда она меняется внутри tree то меняется и z. Такое может быть?

alsav22 · 25.10.2012, 20:40

Сообщение от Default

Такое может быть?

Естественно. В tree() есть вызов funkg() (d = funkg(y,mu,lambda)). После этого память там меняется, а z указатель на эту память. Вам в main() нужно создать массив, и после того, как z получит адрес из funkg(), скопировать значения по этому адресу в массив, чтобы потом вывести эти значения в printf().

Добавлено через 5 минут

C++

double z[2];
...
for(int i=0; i<25;i++)
    {
        x[0] = -0.7+(0.1*i);
        x[1] = -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        double* temp;
        temp = funkg(x, mu, lambda);
        z[0] = temp[0];
        z[1] = temp[1];
        time = clock();
        s4et = tree(x, x0, x1, R0, R1, ispitania1, mu, lambda);
...

@Default · 25.10.2012, 21:04 **[ТС]**

Я так и сделал, но мне кажется, что где-то в других местах тоже может выплыть не то значение, которое нам требуется =(

alsav22 · 26.10.2012, 17:51

Кстати, эта ошибка с z проще устраняется:

C++

int main()
{
    ...
    double *z = NULL;
    ...
 
    for(int i=0; i<25;i++)
    {
        x[0] = -0.7+(0.1*i);
        x[1] = -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        time = clock();
        s4et = tree(x, x0, x1, R0, R1, ispitania1, mu, lambda);
        time = clock() - time;
        z = funkg(x, mu, lambda);
        //fprintf(lame2_1_1,"%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",i+1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
        //if (abs(s4et.answer[0]-z[0])<s4et.dispersion[0]){printf("it is ok");}
        //else {printf("it is not ok");}
        printf("%d %4.1f %8.5f %8.5f %8.5f %8.5f %8.5f %8.5f %d %4.1f\n",1 , x[0], s4et.answer[0], z[0], z[0] - s4et.dispersion[0], z[0] + s4et.dispersion[0], s4et.dispersion[0], s4et.srednee,s4et.max_dlina, (double)time/CLOCKS_PER_SEC);
    }
....

Добавлено через 36 минут
Ещё вариант. Изменить функцию funkg() так, чтобы она заполняла массив, который передаётся ей через указатель. Так как эти массивы используются только внутри функций, то можно массив double g[2] убрать из глобальной области, каждой функции, которая использует вызов funkg(), сделать свой массив и передавать его в funkg() для заполнения. Тогда не будет обращения к общей памяти. Результаты коды выдают идентичные.

C++

...
void funkg(double *g, double koord[],double mu, double lambda);
....
void funkg(double *g, double koord[],double mu, double lambda)//,double Z, double fi)
{
    g[0] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[0],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    g[1] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[1],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
}
...
 
int main()
{
    
    ...
    double z[2];
    ...
    for(int i=0; i<25;i++)
    {
        ...
        funkg(z, x, mu, lambda);
        time = clock();
        s4et = tree(x, x0, x1, R0, R1, ispitania1, mu, lambda);
        time = clock() - time;
        
        ...
        printf("%d %4.1f %8.5f %8.5f %8.5f %8.5f %8.5f %8.5f %d %4.1f\n",1 , x[0], s4et.answer[0], z[0], z[0] - s4et.dispersion[0], z[0] + s4et.dispersion[0], s4et.dispersion[0], s4et.srednee,s4et.max_dlina, (double)time/CLOCKS_PER_SEC);
    }
...
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda)
{
...
double d[2];
...
else
                {
                    funkg(d, y, mu, lambda);
...
else
                {
                    funkg(d, y, mu, lambda);
...

Код

C++

#include <math.h>
#include <time.h>
#include <stdio.h>
#include<iostream>
using namespace std;
 
double pi = 3.14159265358979323846264338327;
long N1rand = 101;
long N2rand = 302;
 
 
double* matr = new double[4]; // указатель, который используется в get_B()
double* del_matr = matr; // указатель для освобождения памяти под matr 
 
double new_matrixS[4]; 
double Ef[4]; 
double matr3[4]; 
 
struct work
{
    double answer[2];
    double dispersion[2];
    double srednee;
    int max_dlina;
};
 
double rand1();
double* direct_multiple(double angle1, double angle2);
void funkg(double *g, double koord[],double mu, double lambda);
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[]);
//изменил имя на matr2, чтобы не совпадало с указателем
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda);
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda);
 
int main()
{
    int ispitania1  = 1000000;
    int ispitania10 = 10000000;
    work s4et;
    double x0[2];
    double x[2];
    double x1[2];
    //double *z = NULL;
    double z[2];
    double R0 = 1;
    double R1 = 1;
    clock_t time;
    double lambda;
    double mu;
    double alpha;
 
    x0[0]=0;
    x0[1]=0;
 
    x1[0]=1;
    x1[1]=1;
 
    mu=1;
    alpha=1;
    lambda=mu*(alpha-1);
    for(int i=0; i<25;i++)
    {
        x[0] = -0.7+(0.1*i);
        x[1] = -0.7+(0.1*i);
        //x[0]= 0.3;//(sqrt(2.0)-1)*cos(pi/4);
        //x[1]= -0.2;//(sqrt(2.0)-1)*sin(pi/4);
        funkg(z, x, mu, lambda);
        time = clock();
        s4et = tree(x, x0, x1, R0, R1, ispitania1, mu, lambda);
        time = clock() - time;
        
        //fprintf(lame2_1_1,"%d %3.1f %3.5f %3.5f %3.5f %3.5f %3.5f %3.5f %d %3.1f\n",i+1, x[0], s4et.answer[0],z[0],z[0]-s4et.dispersion[0],z[0]+s4et.dispersion[0],s4et.dispersion[0],s4et.srednee,s4et.max_dlina,(double)time/CLOCKS_PER_SEC);
        //if (abs(s4et.answer[0]-z[0])<s4et.dispersion[0]){printf("it is ok");}
        //else {printf("it is not ok");}
        printf("%d %4.1f %8.5f %8.5f %8.5f %8.5f %8.5f %8.5f %d %4.1f\n",1 , x[0], s4et.answer[0], z[0], z[0] - s4et.dispersion[0], z[0] + s4et.dispersion[0], s4et.dispersion[0], s4et.srednee,s4et.max_dlina, (double)time/CLOCKS_PER_SEC);
    }
    
    delete [] del_matr;
        
    getchar();
    return 0;
}
 
 
double rand1()/* the initial data(type long):  N1rand=101; N2rand=302; */
{ 
  double x;
  long k;
  k=x=((double)( N1rand=(( N1rand*=1373)%1000919))/1000919.0+
  (double)( N2rand=(( N2rand*=1528)%1400159))/1400159.0 );
  return(x-k);
}
 
 
double* direct_multiple(double angle1, double angle2)
{
    // пока что с двойкой, потом двойка будет входить в переменную для весов
    
    new_matrixS[0]=2*cos(angle1)*cos(angle2);
    new_matrixS[1]=2*cos(angle1)*sin(angle2);
    new_matrixS[2]=2*sin(angle1)*cos(angle2);
    new_matrixS[3]=2*sin(angle1)*sin(angle2);
    return new_matrixS;
}
 
 
void funkg(double *g, double koord[],double mu, double lambda)//,double Z, double fi)
{
    g[0] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[0],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
    g[1] = /*1 + koord[0] - koord[0]*koord[1];//*/1 + pow(koord[1],2) - ((4*mu + 2*lambda)/(mu + lambda))*koord[0]*koord[1];
}
 
 
double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[])
{
    double thetta;
    thetta = atan2(y[1]-Rc[1],y[0]- Rc[0]);
 
    double p0 = 0;
 
    p0 = (alpha*Z - R)/R;
 
    double a1; // случайная величина от 0 до 1 используется для вычисления матрицы B
    double a2; // случайная величина от 0 до 1 используется для вычисления матрицы B если a1 > alpha/(2+alpha)
 
    a1 = rand1();
    a2 = rand1();
 
    if (thetta < 0) thetta = thetta + 2*pi;
    
    if (a1 < alpha/(2+alpha)) matr = direct_multiple(fi,fi);
    else if (2*R - alpha*Z > 0) 
         {
            if (a2 < (2*R - alpha*Z)/2*R)
            {
                Ef[0]=1;
                Ef[1]=0;
                Ef[2]=0;
                Ef[3]=1;
                
                matr = Ef;
            }
            else matr = direct_multiple((thetta + fi)/2,(thetta + fi)/2); 
         }
         else if (a2 < (alpha*Z - 2*R)/(2*R*p0))
              {
                    
                    matr[0] = -1*p0;
                    matr[1] = 0;
                    matr[2] = 0;
                    matr[3] = -1*p0;
              }
              else
              {
                    double* matrthfi = direct_multiple((thetta + fi)/2,(thetta + fi)/2);
                    
                    matr[0] = matrthfi[0]*p0;
                    matr[1] = matrthfi[1]*p0;
                    matr[2] = matrthfi[2]*p0;
                    matr[3] = matrthfi[3]*p0;
              }
    
    return matr;
}
 
 
double* matr2(double t0[],double tn[],double Z,double fi, double radius,double mu, double lambda)
{//t0 - na4alo koord, tn - to4ka reshenia, tr - slu4ainaia to4ka
    
    double costt = (tn[0]-t0[0]+Z*cos(fi))/radius;
 
    double sintt = (tn[1]-t0[1]+Z*sin(fi))/radius;
 
    double cosfi = cos(fi);
 
    double sinfi = sin(fi);
 
    matr3[0]=(((1 + 2*mu/(lambda + mu))-1) + 2*cosfi*cosfi + (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[1]=(2*cosfi*sinfi + (Z/radius)*(sintt*cosfi + costt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    matr3[2]=matr[1];
    
    matr3[3]=(((1 + 2*mu/(lambda + mu))-1) + 2*sinfi*sinfi - (Z/radius)*(costt*cosfi - sintt*sinfi))/(1 + 2*mu/(lambda + mu));
    
    return matr3;
}
 
 
work tree(double tn[],double t0[],double t1[],double R0,double R1,int N, double mu, double lambda)
{
    work resh;
    int k;
    int k2;
    int obr;
    int max_obriv = 0;
    
    double t4[2]; //peremennaia to4ka dlia blujdania
    double y[2];
    int obriv[100];
    
    double E[4];
    double Q[4];
    double Qt[4];
 
    //double* d = NULL;
    double d[2];
    double* matrix = NULL;
 
    double p = 0;
    double modz = 0;
    double N2 = N;
    double fp = 0;
    double srednee = 0;
    int oi = 0;
 
    double p0 = 0; // надо в том случае, если 2R-alpha*Z < 0 
 
    double fi;
    double fipi;
    double J;
    double psi;
    double Z;
    double Zs;
    double a; // alpha случайная величина от 0 до 1
    
    double fis; // угол фи штрих
    double thetta; // угол тетта
 
    double alpha = (mu + lambda)/mu;
 
    resh.answer[0]=0;
    resh.answer[1]=0;
    resh.dispersion[0]=0;
    resh.dispersion[1]=0;
 
    E[0]=1;
    E[1]=0;
    E[2]=0;
    E[3]=1;
 
    for (int pp=0; pp<100; pp++)
        obriv[pp]=0;
 
    if (sqrt(pow((tn[0]-t0[0]),2) + pow((tn[1]-t0[1]),2)) < R0) k2=0;
    else k2=1;
    for(int i=0; i<N; i++)
    {
        obr=0;
        k=k2;
        for(int s=0;s<4;s++)
        {
            Q[s]=E[s];
        }
        for(int s=0;s<2;s++)
        {
            t4[s]=tn[s];
        }
        while(k!=2)
        {
            if (k==0)
            {   
                fis = atan2(t4[1]-t0[1],t4[0]- t0[0]);
                if (fis < 0) {fis = fis + 2*pi;}
                p = sqrt(pow(t4[0]-t0[0],2)+ pow(t4[1]-t0[1],2));
 
                fi = rand1()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R0,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand1();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R0,t0,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
                    
                    matrix = get_B(alpha,fipi,R0,t0,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t1[0]),2) + pow((y[1]-t1[1]),2)) < R1) 
                {
                    k=1;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    //Q=Qt*matr;
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    //continue;
                }
                else
                {
                    funkg(d, y, mu, lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    
                    //break;
                }
                
            } // if (k == 0)
            if (k==1)
            {
                fis = atan2(t4[1]-t1[1],t4[0]- t1[0]);
                p = sqrt(pow(t4[0]-t1[0],2)+ pow(t4[1]-t1[1],2));
 
                fi = rand1()*pi;
                fipi = fi + pi;
                psi = fi - fis;
                J = sqrt(pow(R1,2) - pow(p*sin(psi),2));
                Z = -p*cos(psi) + J;
                Zs = p*cos(psi) + J;
                a = rand1();
 
                if (a <= Zs/(2*J)) // считаем нашу разыграную точку и матрицу переходя для неё первый случай
                {
                    y[0]=t4[0]+Z*cos(fi);
                    y[1]=t4[1]+Z*sin(fi);
 
                    matrix = get_B(alpha,fi,R1,t1,Z,y);
                    
                }
                else // считаем нашу разыграную точку и матрицу переходя для неё второй случай
                {
                    y[0]=t4[0]+Zs*cos(fipi);
                    y[1]=t4[1]+Zs*sin(fipi);
 
                    matrix = get_B(alpha,fipi,R1,t1,Zs,y);
                    
                }
 
                if (sqrt(pow((y[0]-t0[0]),2) + pow((y[1]-t0[1]),2)) < R0)
                {
                    k=0;
                    for(int s=0;s<2;s++)
                    {
                        t4[s]=y[s];
                    }
                    for(int s=0;s<4;s++)
                    {
                        Qt[s]=Q[s];
                    }
                    Q[0]=(Qt[0]*matrix[0]+Qt[1]*matrix[2]);
                    Q[1]=(Qt[0]*matrix[1]+Qt[1]*matrix[3]);
                    Q[2]=(Qt[2]*matrix[0]+Qt[3]*matrix[2]);
                    Q[3]=(Qt[2]*matrix[1]+Qt[3]*matrix[3]);
 
                    obr++;
                    
                    //continue;
                }
                else
                {
                    funkg(d, y, mu, lambda);
                    resh.answer[0]+=((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]);
                    resh.dispersion[0]+=pow(((Q[0]*matrix[0]+Q[1]*matrix[2])*d[0]+(Q[0]*matrix[1]+Q[1]*matrix[3])*d[1]),2);
                    resh.answer[1]+=((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]);
                    resh.dispersion[1]+=pow(((Q[2]*matrix[0]+Q[3]*matrix[2])*d[0]+(Q[2]*matrix[1]+Q[3]*matrix[3])*d[1]),2);
                    k=2;
                    //break;
                }
            } //if (k ==1)
                
        } // while
    if (obr < 100) obriv[obr]++;
    }
    
    int obr2 = 99;
    while(obriv[obr2]==0) obr2--;
    resh.max_dlina = obr2+1;    
 
    for (int i=0; i<100; i++)
    {
        srednee+=(i+1)*obriv[i];
    }
    resh.srednee = srednee/N2;
        
    resh.answer[0]=resh.answer[0]/N;
    resh.answer[1]=resh.answer[1]/N;
    
    resh.dispersion[0]= 3*sqrt(resh.dispersion[0]/N - pow(resh.answer[0],2)) / sqrt(N2);
    resh.dispersion[1]= 3*sqrt(resh.dispersion[1]/N - pow(resh.answer[1],2)) / sqrt(N2);
 
    return resh;
}

Добавлено через 9 минут
Ещё вот Ef[4]. Используется в одном месте (в get_B()) и всегда с одинаковыми значениями. Не проще ли тогда убрать этот массив, а в get_B() так сделать:

C++

double* get_B(double alpha,double fi, double R, double Rc[], double Z, double y[])
{
    ...
            if (a2 < (2*R - alpha*Z)/2*R)
            {
                /*Ef[0]=1;
                Ef[1]=0;
                Ef[2]=0;
                Ef[3]=1;*/
 
                matr[0]=1;
                matr[1]=0;
                matr[2]=0;
                matr[3]=1;
                
                //matr = Ef;
            }
            else matr = direct_multiple((thetta + fi)/2,(thetta + fi)/2); 
        ...
    return matr;
}

И советую изменить имя у функции rand(). Есть стандартная функция с таким именем, и если подключать другие библиотеки, то возникает конфликт имён.

Croessmah Неэпический 17869 / 10634 / 2054 Регистрация: 27.09.2012 Сообщений: 26,736 Записей в блоге: 1
	24.10.2012, 18:41	21
	Сообщение от Default Как я догадываюсь основная проблема у меня в том, что когда я из функции get_B() вызываю функцию direct_multiple() у меня происходит накопление памяти. Как от него избавиться не ясно. Когда перед return прописываю delete[] matrfi и matrthfi, в которые я передаю direct_multiple() вылазит фатал, из-за чего конкретно и в чем косяк разобраться пока что не выходит. Исходя из сказанного Вами, могу сделать вывод - Вы просто не понимаете как работать с памятью. Посему еще раз говорю - откройте книжку и почитайте. Можете взять любую книгу по C++ и почитать главу "Динамическая память(куча)" 0

@Default 0 / 0 / 0 Регистрация: 08.04.2012 Сообщений: 28
	24.10.2012, 18:43 [ТС]	22
	Не понимаю, сейчас буду смотреть видео лекцию / искать книги. 0

@Default 0 / 0 / 0 Регистрация: 08.04.2012 Сообщений: 28
	25.10.2012, 09:25 [ТС]	24
	Огромное спасибо, теперь буду вникать и пытаться разбираться как оно теперь работает =) 0

@Default 0 / 0 / 0 Регистрация: 08.04.2012 Сообщений: 28
	25.10.2012, 19:45 [ТС]	26
	Появилась проблема, которую я сейчас пытаюсь отловить, в мэйне z=funkg(x,mu,lambda); time = clock(); s4et = tree(x,x0,x1,R0,R1,ispitania1,mu,lambda); time = clock() - time; после вызова функции tree, меняется значение z, что никак не допустимо для вычислений. Видимо z идет как указатель на funkg() и когда она меняется внутри tree то меняется и z. Такое может быть? 0

@Default 0 / 0 / 0 Регистрация: 08.04.2012 Сообщений: 28
	25.10.2012, 21:04 [ТС]	28
	Я так и сделал, но мне кажется, что где-то в других местах тоже может выплыть не то значение, которое нам требуется =( 0

	26.10.2012, 17:51