Форум программистов, компьютерный форум CyberForum.ru

Коллекция алгоритмов от Johna Smith - C++

Восстановить пароль Регистрация
 
Рейтинг: Рейтинг темы: голосов - 17, средняя оценка - 5.00
LK
Заблокирован
19.05.2013, 14:32     Коллекция алгоритмов от Johna Smith #1
Коллекция алгоритмов от Johna Smith
в качестве учебного материала

Источник: http://vingrad.ru
Выложенные здесь алгоритмы преследуют исключительно учебные цели.
Код неоптимизирован, местами морально устарел и показывает только принцип решения той или иной задачи. Однако он вполне рабочий и может быть использован с соответствующей дорботкой под Ваши собственные нужды.
схема, не для копипейста, не проверял
консультации и техническая поддержка не предоставляются


Добавлено через 4 минуты
Перевод из... в...

Градусы, минуты, секунды - в градусы
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
//////////////////////////////////////////////////////////////////////////////
//
//  Translating angle from degrees, minutes, seconds to degrees
//  (c) Johna Smith, 1996
//
//  Method description:
//    Here we just use the following formula:
//      phi=(phi"/60+phi')/60+phiш
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
struct angle
{
  float degrees;
  float minutes;
  float seconds;
};
void print_angle(angle a)
{
  printf("%fш%f'%f\"",a.degrees,a.minutes,a.seconds);
}
float dms_to_d(angle a)
{
  float f;
  f=(a.seconds/60+a.minutes)/60+a.degrees;
  return f;
}
void main(void)
{
  angle a={30,30,30};
  print_angle(a);
  printf("= %fш",dms_to_d(a));
}

Градусы - в градусы, минуты, секунды
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
//////////////////////////////////////////////////////////////////////////////
//
//  Translating angle from degrees to degrees, minutes, seconds
//  (c) Johna Smith, 1996
//
//  Method description:
//    degrees is the integer part of angle f
//    minutes is the integer part of the remainder multiplied by 60
//    seconds is the integer part of 60*(phi_in_minutes-phi')
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
struct angle
{
  float degrees;
  float minutes;
  float seconds;
};
void print_angle(angle a)
{
  printf("%fш%f'%f\"",a.degrees,a.minutes,a.seconds);
}
angle d_to_dms(float f)
{
  angle a;
  a.degrees=(int)f;
  a.minutes=(int)((f-a.degrees)*60);
  a.seconds=(int)((f-a.degrees)*60-a.minutes);
  return a;
}
void main(void)
{
  printf("12.3456ш=");
  print_angle(d_to_dms(12.3456));
}

Комплексных величин - в экспоненциальную форму
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (converting from a+bi form to M*exp(i*phi))
//  (c) Johna Smith, 1996
//
//  Given: z - complex value
//  z=a+i*b
//  z=M*exp(i*phi)
//  M=(a^2+b^2)^(1/2)    phi=arccos(a/|M|)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
struct exp_complex
{
  float M;
  float phi;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
void show_exp_complex(exp_complex c) // this functions displays complex value
{
  printf("%f*exp(%f*i)",c.M,c.phi);
}
exp_complex Convert(complex a)
{
  exp_complex b;
  b.M=sqrt(a.re*a.re+a.im*a.im);
  b.phi=(a.im<0?-1:1)*acos(a.re/fabs(b.M));
  return b;
}
complex a={-1,1};
void main(void)
{
  show_complex(a);
  printf(" = ");
  show_exp_complex(Convert(a));
}

Градусы по Фаренгейту - в градусы по Цельсию
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
//////////////////////////////////////////////////////////////////////////////
//
//  Translatintg Farengheit degrees to Celsy degrees and in other direction
//  (c) Johna Smith, 1996
//
//  Method description:
//
//    oF = 5/9(n-32) oC     oC = (32+9n/5) oF
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
float Cels2Fareng(float degree)
{
  return (5*(degree-32)/9);
}
float Fareng2Cels(float degree)
{
  return (32+9*degree/5);
}
void main(void)
{
  printf("100oC = %f oF\n",Cels2Fareng(100));
  printf("-50oF = %f oC\n",Fareng2Cels(-50));
}

Декартовы координаты - в полярные (2D)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Translatintg Decart coordinates to polar coordinates in 2 dimensions
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  r=(x^2+y^2)^(1/2)
//  phi=+/- arccos(x/r)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct decart
{
  float x;
  float y;
};
struct polar
{
  float r;
  float phi;
};
polar Decart2Polar(decart c)
{
  polar p;
  p.r=sqrt(c.x*c.x+c.y*c.y);
  p.phi=acos(c.x/p.r)*(c.y>=0?1:-1);
  return p;
}
decart d={1,1};
polar p;
void main(void)
{
  p=Decart2Polar(d);
  printf("(x,y)=(%f,%f) -> (r,phi)=(%f,%f)\n",d.x,d.y,p.r,p.phi);
}

Полярные - в декартовы
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
//////////////////////////////////////////////////////////////////////////////
//
//  Translatintg polar coordinates to Decart coordinates in 2 dimensions
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  x=r*cos(phi)
//  y=r*sin(phi)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi 3.1415926536
struct decart
{
  float x;
  float y;
};
struct polar
{
  float r;
  float phi;
};
decart Polar2Decart(polar p)
{
  decart d;
  d.x=p.r*cos(p.phi);
  d.y=p.r*sin(p.phi);
  return d;
}
polar p={1.4142135624,Pi/4};
decart d;
void main(void)
{
  d=Polar2Decart(p);
  printf("(r,phi)=(%f,%f) -> (x,y)=(%f,%f)\n",p.r,p.phi,d.x,d.y);
}

Декартовы координаты - в сферические (3D)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
//////////////////////////////////////////////////////////////////////////////
//
//  Translatintg Decart coordinates to spherical coordinates in 3 dimensions
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  r=(x^2+y^2+z^2)^(1/2)
//  phi=+/- arccos(x/(x*x+y*y)^(1/2))
//  theta=arccos(z/r)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct decart
{
  float x;
  float y;
  float z;
};
struct spherical
{
  float r;
  float phi;
  float theta;
};
spherical Decart2Spherical(decart c)
{
  spherical p;
  p.r=sqrt(c.x*c.x+c.y*c.y+c.z*c.z);
  p.phi=acos(c.x/sqrt(c.x*c.x+c.y*c.y))*(c.y>=0?1:-1);
  p.theta=acos(c.z/p.r);
  return p;
}
decart d={1,1,1};
spherical p;
void main(void)
{
  p=Decart2Spherical(d);
  printf("(x,y,z)=(%f,%f,%f) -> (r,phi,theta)=(%f,%f,%f)\n",
         d.x,d.y,d.z,p.r,p.phi,p.theta);
}

Сферические - в декартовы
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
//////////////////////////////////////////////////////////////////////////////
//
//  Translatintg spherical coordinates to Decart coordinates in 3 dimensions
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  x=r*cos(phi)*sin(theta)
//  y=r*sin(phi)*sin(theta)
//  z=r*cos(theta)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi 3.1415926536
struct decart
{
  float x;
  float y;
  float z;
};
struct spherical
{
  float r;
  float phi;
  float theta;
};
decart Spherical2Decart(spherical p)
{
  decart d;
  d.x=p.r*cos(p.phi)*sin(p.theta);
  d.y=p.r*sin(p.phi)*sin(p.theta);
  d.z=p.r*cos(p.theta);
  return d;
}
spherical p={1.732051,0.785398,0.955317};
decart d;
void main(void)
{
  d=Spherical2Decart(p);
  printf("(r,phi,theta)=(%f,%f,%f) -> (x,y,z)=(%f,%f,%f)\n",
         p.r,p.phi,p.theta,d.x,d.y,d.z);
}

Число из системы с основанием М - в систему с основанием N
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
//////////////////////////////////////////////////////////////////////////////
//
//  Translating number from M-system to N-system
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) translate number to decimal
//    2) translate from decimal using sequence of divisions by nbase
//       remainders will be digits of the number in the new system
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
enum digit {A=10,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z};
char mbase=16;    // from hex
char nbase=2;  // to binary
digit m[]={F,E,D,C,B,A};  // (FEDCBA)(16)
digit n[255];
void main(void)
{
  unsigned long int decimal=0,tmp=1;
  float tmpf;
  // printing given number
  for (int i=0;i<sizeof(m)/sizeof(digit);i++)
    if (m[i]<10)
       printf("%d",m[i]);
    else printf("%c",m[i]+55);
  
  // translate number to decimal
  for(i=sizeof(m)/sizeof(digit)-1;i>=0;i--)
  {
     decimal+=tmp*m[i];
     tmp*=mbase;
  }
  printf("(%d) = %ld(10) = ",mbase,decimal);
  
  // translating from decimal
  tmpf=(float)log(decimal/(nbase-1))/log(nbase);
  // rounding to nearest integer
  tmp=(int)tmpf;
  
  if (tmpf-tmp>0.5) tmp++;
  
  for (i=tmp;i>=0;i--)
  {
    n[i]=decimal%nbase;
    decimal/=nbase;
  }
  
  // printing converted number
  for (i=0;i<=tmp;i++)
    if (n[i]<10)
      printf("%d",n[i]);
    else printf("%c",n[i]+55);
  printf("(%d)\n",nbase);
}

Приведение периодической дробной части к нормальному виду
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
//////////////////////////////////////////////////////////////////////////////
//
//  Transforming periodical decimal fraction to normal one.
//  (c) Johna Smith, 1996
//
//  Method description:
//    This method uses the fact that periodical fraction can be
//   calculated as a sum of infinite geometr. progression = 1/(1-q).
//    For example: 235.353535353...
//    The period is 35, the base (aperiodical part of the fraction) is 200,
//   number of digits in the period is 2. These parameters are required by
//   function Convert.
//    q=10^(-2)  x0=35   SUM=x0+q*x0+q^2*x0+q^3*x0+...
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct fraction
{
  long int a;
  long int b;
  long int c;
} f;    // f=a+b/c
fraction Convert(double base, double period, int digits)
{
  fraction f;
  double tmp;
  f.a=0;
  f.b=pow10(digits);
  f.c=f.b-1;
  // calculating 1/(1-q)
  tmp=(double)f.b*period;
  while (floor(tmp)!=tmp)
  {
     tmp*=10;
     f.b*=10;
     f.c*=10;
  }
  f.b*=period;
  // adding base to the fraction
  base*=f.c;
  while (floor(base)!=base)
  {
     base*=10;
     f.b*=10;
     f.c*=10;
  }
  f.b+=base;
  // if this fraction isn't right: b>c we should convert it
  // into the right fraction a+ b1/c, where b1<c1 & a*c+b1=b
  if (f.b/f.c>1)
  {
     f.a+=(int)(f.b/f.c);
     f.b-=((int)(f.b/f.c))*f.c;
  }
  return f;
}
void main(void)
{
  f=Convert(200,35,2);
  printf("235,(35) = %ld %ld/%ld\n",f.a,f.b,f.c);
  f=Convert(0.0783,441e-7,3);
  printf("0.0783(441) = %ld %ld/%ld",f.a,f.b,f.c);
}


Добавлено через 18 минут
Графические алгоритмы

Рисование 3D объекта
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
//////////////////////////////////////////////////////////////////////////////
//
//  Drawing 3D object
//  (c) Johna Smith, 1996
//
//  Method description:
//    Function draw3Dobject recieves the following parameters:
//     object - reference to array of points of the drawn object
//     N - number of points in this array
//     rho     \
//     phi     |- spherical coordinates of point of view
//     theta   /
//     dist_to_screen - distance from viewpoint to screen
//     xshift, yshift - picture will be shifted by this values
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
#define Pi 3.1415926536
enum Action{move,draw};
struct Point3D
{
  int x;
  int y;
  int z;
  Action action;
};
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
        /* request autodetection */
        int gdriver = DETECT, gmode, errorcode;
        /* initialize graphics mode */
        initgraph(&gdriver, &gmode, "");
        /* read result of initialization */
        errorcode = graphresult();
        if (errorcode != grOk)    /* an error occurred */
      {
                printf("Graphics error: %s\n", grapherrormsg(errorcode));
                printf("Press any key to halt:");
                                         getch();
                exit(1);               /* return with error code */
        }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
void draw3Dobject(Point3D *object, int N, float rho, float theta,
                  float phi, float dist_to_screen, int xshift, int yshift)
{
  int x,y;
  float xe,ye,ze,costh,sinph,cosph,sinth,v11,v12,v13,v21,v22,v32,v33,v23,v43;
  // calculating coefficients
  costh=cos(theta);
  sinth=sin(theta);
  cosph=cos(phi);
  sinph=sin(phi);
  v11=-sinth; v12=-cosph*costh; v13=-sinph*costh;
  v21=costh;  v22=-cosph*sinth; v23=-sinph*sinth;
                                  v32=sinph;        v33=-cosph;
                                v43=rho;
  for (int i=0;i<N;i++)
  {
      // calculating eye coordinates
      xe=v11*(object+i)->x+v21*(object+i)->y;
      ye=v12*(object+i)->x+v22*(object+i)->y+v32*(object+i)->z;
      ze=v13*(object+i)->x+v23*(object+i)->y+v33*(object+i)->z+v43;
      // calculating screen coordinates
      x=dist_to_screen*xe/ze+xshift;
      y=dist_to_screen*ye/ze+yshift;
      // drawing
      if((object+i)->action==move)
          MoveTo(x,y);
      else
          LineTo(x,y);
  }
}
int main(void)
{
  Point3D thetr[]={{100,100,100,move},
                              {100,200,100,draw},
                              {200,100,100,draw},
                              {100,100,100,draw},
                              {100,100,200,draw},
                              {200,100,100,draw},
                              {100,100,200,move},
                              {100,200,100,draw}};
  int N=sizeof(thetr)/sizeof(Point3D);
  float rho=700,theta=Pi/4,phi=Pi/4,dist_to_screen=300;
  int xshift=300, yshift=150;
  // initializing graphics mode
  init_gr();
  /* examples */
  while (!kbhit())
  {
      theta+=0.05;
      // rotating viewpoint
      if (phi>2*Pi) phi-=2*Pi;
      setcolor(11);
      draw3Dobject(thetr,N,rho,theta,phi,dist_to_screen,xshift,yshift);
      setcolor(0);
      delay(15);
      draw3Dobject(thetr,N,rho,theta,phi,dist_to_screen,xshift,yshift);
  }
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Вращение 3D объекта
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
//////////////////////////////////////////////////////////////////////////////
//
//  Rotating 3D object
//  (c) Johna Smith, 1996
//
//  Method description:
//   Given: vector X, angle a
//   To rotate this vector about axe X we should apply the following
//   operation:  X'=X*Rx, where Rx - is matrix of rotation
//        [ 1 0      0     0 ]
//        [ 0 cos a  sin a 0 ]
//     Rx=[ 0 -sin a cos a 0 ]
//        [ 0 0      0     1 ]
//   Applying this operation to every point of the given object we
//   rotate it.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
#define Pi 3.1415926536
enum Action{move,draw};
struct Point3D
{
  float x;
  float y;
  float z;
  Action action;
};
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
        /* request autodetection */
        int gdriver = DETECT, gmode, errorcode;
        /* initialize graphics mode */
        initgraph(&gdriver, &gmode, "");
        /* read result of initialization */
        errorcode = graphresult();
        if (errorcode != grOk)    /* an error occurred */
        {
               printf("Graphics error: %s\n", grapherrormsg(errorcode));
               printf("Press any key to halt:");
               getch();
               exit(1);               /* return with error code */
        }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
void draw3Dobject(Point3D *object, int N, float rho, float theta,
                  float phi, float dist_to_screen, int xshift, int yshift)
{
  int x,y;
  float xe,ye,ze,costh,sinph,cosph,sinth,v11,v12,v13,v21,v22,v32,v33,v23,v43;
  // calculating coefficients
  costh=cos(theta);
  sinth=sin(theta);
  cosph=cos(phi);
  sinph=sin(phi);
  v11=-sinth; v12=-cosph*costh; v13=-sinph*costh;
  v21=costh;  v22=-cosph*sinth; v23=-sinph*sinth;
                      v32=sinph;             v33=-cosph;
                                                     v43=rho;
  for (int i=0;i<N;i++)
  {
         // calculating eye coordinates
         xe=v11*(object+i)->x+v21*(object+i)->y;
         ye=v12*(object+i)->x+v22*(object+i)->y+v32*(object+i)->z;
         ze=v13*(object+i)->x+v23*(object+i)->y+v33*(object+i)->z+v43;
         // calculating screen coordinates
         x=dist_to_screen*xe/ze+xshift;
         y=dist_to_screen*ye/ze+yshift;
 
        // drawing
         if((object+i)->action==move)
              MoveTo(x,y);
         else
              LineTo(x,y);
  }
}
void RotateX(Point3D *object, int n, float angle)
{
  float cosa,sina,y,z;
  cosa=cos(angle);
  sina=sin(angle);
  for(int i=0;i<n;i++)
  {
         y=(object+i)->y*cosa-(object+i)->z*sina;
         z=(object+i)->y*sina+(object+i)->z*cosa;
         (object+i)->y=y;
         (object+i)->z=z;
  }
}
void RotateY(Point3D *object, int n, float angle)
{
  float cosa,sina,x,z;
  cosa=cos(angle);
  sina=sin(angle);
  for(int i=0;i<n;i++)
  {
         x=(object+i)->x*cosa+(object+i)->z*sina;
         z=-(object+i)->x*sina+(object+i)->z*cosa;
         (object+i)->x=x;
         (object+i)->z=z;
  }
}
void RotateZ(Point3D *object, int n, float angle)
{
  float cosa,sina,x,y;
  cosa=cos(angle);
  sina=sin(angle);
  for(int i=0;i<n;i++)
  {
         x=(object+i)->x*cosa-(object+i)->y*sina;
         y=(object+i)->x*sina+(object+i)->y*cosa;
         (object+i)->x=x;
         (object+i)->y=y;
  }
}
int main(void)
{
  Point3D thetr[]={{100,100,100,move},
                              {100,200,100,draw},
                              {200,100,100,draw},
                              {100,100,100,draw},
                              {100,100,200,draw},
                              {200,100,100,draw},
                              {100,100,200,move},
                              {100,200,100,draw}};
  Point3D cube[]={{-50,-50,-50,move},
                              {-50,50,-50,draw},
                              {-50,50,50,draw},
                              {50,50,50,draw},
                              {50,50,-50,draw},
                              {50,-50,-50,draw},
                              {50,-50,50,draw},
                              {-50,-50,50,draw},
                              {-50,50,50,draw},
                              {-50,-50,50,move},
                              {-50,-50,-50,draw},
                              {50,-50,-50,draw},
                              {50,50,-50,move},
                              {-50,50,-50,draw},
                              {50,-50,50,move},
                              {50,50,50,draw}};
  int N=sizeof(thetr)/sizeof(Point3D);
  int M=sizeof(cube)/sizeof(Point3D);
  float rho=1500,theta=Pi/4,phi=-3*Pi/4,dist_to_screen=700;
  int xshift=300, yshift=150,xshift1=200,yshift1=200;
  // initializing graphics mode
  init_gr();
  /* examples */
  while (!kbhit())
  {
          RotateX(cube,M,Pi/24);
          RotateY(cube,M,Pi/24);
          RotateZ(cube,M,Pi/24);
          RotateY(thetr,N,Pi/30);
          setcolor(12);
          draw3Dobject(thetr,N,rho,theta,phi,dist_to_screen,xshift,yshift);
          draw3Dobject(cube,M,rho,theta,phi,dist_to_screen,xshift1,yshift1);
          setcolor(0);
          delay(35);
          draw3Dobject(thetr,N,rho,theta,phi,dist_to_screen,xshift,yshift);
          draw3Dobject(cube,M,rho,theta,phi,dist_to_screen,xshift1,yshift1);
  }
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Рисование куба
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
//////////////////////////////////////////////////////////////////////////////
//
//  Generating and drawing cube
//  (c) Johna Smith, 1996
//
//  Method description:
//    Cube is one of simplest figures in stereometry.
//    We just need to generate 'n' squares parallel to X,Y and Z axes
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
#define Pi 3.1415926536
enum Action{move,draw};
struct Point3D
{
  int x;
  int y;
  int z;
  Action action;
};
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
  /* request autodetection */
  int gdriver = DETECT, gmode, errorcode;
  /* initialize graphics mode */
  initgraph(&gdriver, &gmode, "");
  /* read result of initialization */
  errorcode = graphresult();
  if (errorcode != grOk)    /* an error occurred */
  {
     printf("Graphics error: %s\n", grapherrormsg(errorcode));
     printf("Press any key to halt:");
     getch();
     exit(1);               /* return with error code */
  }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
void draw3Dobject(Point3D *object, int N, float rho, float theta,
                  float phi, float dist_to_screen, int xshift, int yshift)
{
  int x,y;
  float xe,ye,ze,costh,sinph,cosph,sinth,v11,v12,v13,v21,v22,v32,v33,v23,v43;
  // calculating coefficients
  costh=cos(theta);
  sinth=sin(theta);
  cosph=cos(phi);
  sinph=sin(phi);
  v11=-sinth; v12=-cosph*costh; v13=-sinph*costh;
  v21=costh;  v22=-cosph*sinth; v23=-sinph*sinth;
              v32=sinph;        v33=-cosph;
                                v43=rho;
  for (int i=0;i<N;i++)
  {
     // calculating eye coordinates
     xe=v11*(object+i)->x+v21*(object+i)->y;
     ye=v12*(object+i)->x+v22*(object+i)->y+v32*(object+i)->z;
     ze=v13*(object+i)->x+v23*(object+i)->y+v33*(object+i)->z+v43;
   
     // calculating screen coordinates
     x=dist_to_screen*xe/ze+xshift;
     y=dist_to_screen*ye/ze+yshift;
     // drawing
     if((object+i)->action==move)
       MoveTo(x,y);
     else
       LineTo(x,y);
  }
}
int main(void)
{
  const int n=4; // number of cubes segments +1
  Point3D cube[15*n]; // coordinates for cubes points
  float rho=1900,theta=Pi/3,phi=2*Pi/3,dist_to_screen=600; // view point
  int xshift=300, yshift=140; // picture offset
  float side=300; // cubes parameters
  float delta;// auxulary variable
  // initializing graphics mode
  init_gr();
  // generating cube
  delta=side/(n-1);
  for (int i=0;i<n;i++)
  {
    for(int j=i*5;j<i*5+5;j++)
    {
      cube[j].x=i*delta;
      cube[j].action=draw;
    }
    cube[i*5].y=0;
    cube[i*5].z=0;
    cube[i*5].action=move;
    cube[i*5+1].y=side;
    cube[i*5+1].z=0;
    cube[i*5+2].y=side;
    cube[i*5+2].z=side;
    cube[i*5+3].y=0;
    cube[i*5+3].z=side;
    cube[i*5+4].y=0;
    cube[i*5+4].z=0;
  }
  int c=5*n;
  for (i=0;i<n;i++)
  {
    for(int j=i*5;j<i*5+5;j++)
    {
      cube[c+j].y=i*delta;
      cube[c+j].action=draw;
    }
    cube[c+i*5].x=0;
    cube[c+i*5].z=0;
    cube[c+i*5].action=move;
    cube[c+i*5+1].x=side;
    cube[c+i*5+1].z=0;
    cube[c+i*5+2].x=side;
    cube[c+i*5+2].z=side;
    cube[c+i*5+3].x=0;
    cube[c+i*5+3].z=side;
    cube[c+i*5+4].x=0;
    cube[c+i*5+4].z=0;
  }
  c=10*n;
  for (i=0;i<n;i++)
  {
    for(int j=i*5;j<i*5+5;j++)
    {
       cube[c+j].z=i*delta;
       cube[c+j].action=draw;
    }
    cube[c+i*5].y=0;
    cube[c+i*5].x=0;
    cube[c+i*5].action=move;
    cube[c+i*5+1].y=side;
    cube[c+i*5+1].x=0;
    cube[c+i*5+2].y=side;
    cube[c+i*5+2].x=side;
    cube[c+i*5+3].y=0;
    cube[c+i*5+3].x=side;
    cube[c+i*5+4].y=0;
    cube[c+i*5+4].x=0;
  }
  // drawing
  draw3Dobject(cube,15*n,rho,theta,phi,dist_to_screen,xshift,yshift);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Рисование тора
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
//////////////////////////////////////////////////////////////////////////////
//
//  Generating and drawing torus
//  (c) Johna Smith, 1996
//
//  Method description:
//    Torus has two main circles, which are described by
//   two systems of eqations:
//    x=Rcos(a)     x=R+rcos(b)
//    y=Rsin(a)     y=0
//    z=0           z=rsin(b)
//    By generating this circles (approximating by polygons) we can get torus
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
#define Pi 3.1415926536
enum Action{move,draw};
struct Point3D
{
  int x;
  int y;
  int z;
  Action action;
};
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
void draw3Dobject(Point3D *object, int N, float rho, float theta,
                  float phi, float dist_to_screen, int xshift, int yshift)
{
  int x,y;
  float xe,ye,ze,costh,sinph,cosph,sinth,v11,v12,v13,v21,v22,v32,v33,v23,v43;
  // calculating coefficients
  costh=cos(theta);
  sinth=sin(theta);
  cosph=cos(phi);
  sinph=sin(phi);
  v11=-sinth; v12=-cosph*costh; v13=-sinph*costh;
  v21=costh;  v22=-cosph*sinth; v23=-sinph*sinth;
              v32=sinph;        v33=-cosph;
                                v43=rho;
  for (int i=0;i<N;i++)
  {
      // calculating eye coordinates
      xe=v11*(object+i)->x+v21*(object+i)->y;
      ye=v12*(object+i)->x+v22*(object+i)->y+v32*(object+i)->z;
      ze=v13*(object+i)->x+v23*(object+i)->y+v33*(object+i)->z+v43;
 
      // calculating screen coordinates
      x=dist_to_screen*xe/ze+xshift;
      y=dist_to_screen*ye/ze+yshift;
   
      // drawing
      if((object+i)->action==move)
        MoveTo(x,y);
      else
        LineTo(x,y);
  }
}
int main(void)
{
  const int n=20;  // number of torus' segments
  Point3D torus[2*n*(n+1)]; // coordinates for torus' points
  float rho=1800,theta=0,phi=3*Pi/4,dist_to_screen=600; // view point
  int xshift=300, yshift=250; // picture offset
  float delta=2.0*Pi/n, r=75, R=300; // torus' parameters
  float alpha,cosa,sina,beta,x; // auxulary variables
  // initializing graphics mode
  init_gr();
  // generating torus
  for (int i=0;i<n;i++)
  {
     alpha=i*delta;
     cosa=cos(alpha);
     sina=sin(alpha);
     for (int j=0;j<n+1;j++)
     {
       beta=j*delta;
       x=R+r*cos(beta);
       torus[i*(n+1)+j].x=cosa*x;
       torus[i*(n+1)+j].y=sina*x;
       torus[i*(n+1)+j].z=r*sin(beta);
       torus[i*(n+1)+j].action=((i==0 && j==0)?move:draw);
     }
  }
  int c=n*n+n;
  for (i=0;i<n;i++)
  {
     beta=i*delta;
     x=R+r*cos(beta);
     for (int j=0;j<n+1;j++)
     {
       alpha=j*delta;
       cosa=cos(alpha);
       sina=sin(alpha);
       torus[c+i*(n+1)+j].x=cosa*x;
       torus[c+i*(n+1)+j].y=sina*x;
       torus[c+i*(n+1)+j].z=r*sin(beta);
       torus[c+i*(n+1)+j].action=draw;
     }
  }
  // drawing
  draw3Dobject(torus,2*n*(n+1),rho,theta,phi,dist_to_screen,xshift,yshift);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Рисование полусферы
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
//////////////////////////////////////////////////////////////////////////////
//
//  Generating and drawing semisphere
//  (c) Johna Smith, 1996
//
//  Method description:
//    Hemisphere is described by following eqation:
//       2    2    2
//      x  + y  + z  = 1,    -1<=z<=0
//
//    By generating this points with step 'delta' we can approximate hemisphere
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
#define Pi 3.1415926536
enum Action{move,draw};
struct Point3D
{
  int x;
  int y;
  int z;
  Action action;
};
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
void draw3Dobject(Point3D *object, int N, float rho, float theta,
                    float phi, float dist_to_screen, int xshift, int yshift)
{
  int x,y;
  float xe,ye,ze,costh,sinph,cosph,sinth,v11,v12,v13,v21,v22,v32,v33,v23,v43;
  // calculating coefficients
  costh=cos(theta);
  sinth=sin(theta);
  cosph=cos(phi);
  sinph=sin(phi);
  v11=-sinth; v12=-cosph*costh; v13=-sinph*costh;
  v21=costh;  v22=-cosph*sinth; v23=-sinph*sinth;
              v32=sinph;        v33=-cosph;
                                v43=rho;
  for (int i=0;i<N;i++)
  {
     // calculating eye coordinates
     xe=v11*(object+i)->x+v21*(object+i)->y;
     ye=v12*(object+i)->x+v22*(object+i)->y+v32*(object+i)->z;
     ze=v13*(object+i)->x+v23*(object+i)->y+v33*(object+i)->z+v43;
    
     // calculating screen coordinates
     x=dist_to_screen*xe/ze+xshift;
     y=dist_to_screen*ye/ze+yshift;
     // drawing
     if((object+i)->action==move)
       MoveTo(x,y);
     else
       LineTo(x,y);
  }
}
int main(void)
{
  const int n=10;  // number of hemispheres segments
  Point3D hemisphere[8*n*n+n]; // coordinates for hemispheres points
  float rho=1800,theta=Pi,phi=3*Pi/4,dist_to_screen=600; // view point
  int xshift=300, yshift=250; // picture offset
  float delta=Pi/(2*n), R=300; // hemisphere parameters
  float alpha,cosa,sina,beta,cosb,sinb; // auxulary variables
  // initializing graphics mode
  init_gr();
  // generating semisphere
  for (int i=0;i<4*n;i++)
  {
     alpha=i*delta;
     cosa=cos(alpha);
     sina=sin(alpha);
     for (int j=0;j<n;j++)
     {
       beta=j*delta;
       sinb=sin(beta);
       cosb=cos(beta);
       hemisphere[i*n+j].x=R*cosa*sinb;
       hemisphere[i*n+j].y=R*sina*sinb;
       hemisphere[i*n+j].z=-R*cosb;
       hemisphere[i*n+j].action=(j==0?move:draw);
     }
  }
  int c=4*n*n;
  for (i=0;i<n;i++)
  {
     beta=i*delta;
     sinb=sin(beta);
     cosb=cos(beta);
     for (int j=0;j<4*n+1;j++)
     {
       alpha=j*delta;
       cosa=cos(alpha);
       sina=sin(alpha);
       hemisphere[c+i*(4*n+1)+j].x=R*cosa*sinb;
       hemisphere[c+i*(4*n+1)+j].y=R*sina*sinb;
       hemisphere[c+i*(4*n+1)+j].z=-R*cosb;
       hemisphere[c+i*(4*n+1)+j].action=(j==0?move:draw);
     }
  }
  // drawing
  draw3Dobject(hemisphere,8*n*n+n,rho,theta,phi,dist_to_screen,xshift,yshift);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Вращение фигуры в плоскости
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
//////////////////////////////////////////////////////////////////////////////
//
//  Rotations in 2 dimensions
//  (c) Johna Smith, 1996
//
//  Method description:
//    There is the following formulas for rotating in 2 dimensions:
//      x'=x0+(x-x0)*cos(phi)-(y-y0)*sin(phi)
//      y'=y0+(x-x0)*sin(phi)+(y-y0)*cos(phi)
//    To rotate the figure we should rotate each of its points
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#define Pi 3.1415926536
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
const N=6; // number of points in the figure
// coordinates of all given points
enum actions {MOVE,DRAW};
struct
{
  actions action;
  int x;
  int y;
} figure[N]={{MOVE,360,270},{DRAW,360,260},{DRAW,355,260},{DRAW,360,250},
     {DRAW,365,260},{DRAW,360,260}};
int x0,y0,dx,dy;
float phi;
int main(void)
{
  // initializing graphics mode
  init_gr();
  // rotating about (x0,y0)
  x0=300;
  y0=260;
  // by 10 degrees
  phi=45.0*Pi/180.0;
  // main loop
  for(int i=0;i<8;i++)
  {
    // rotating the figure
    for (int j=0;j<N;j++)
    {
      dx=figure[j].x-x0;
      dy=figure[j].y-y0;
      figure[j].x=x0+dx*cos(phi)-dy*sin(phi);
      figure[j].y=y0+dx*sin(phi)+dy*cos(phi);
    }
    // drawing rotated figure
  for (j=0;j<N;j++)
      if (figure[j].action==MOVE)
        MoveTo(figure[j].x,figure[j].y);
      else
        LineTo(figure[j].x,figure[j].y);
  }
  // clean up
  getch();
  end_gr();
  return 0;
}

Рисование линии (по Брезенхэму)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
//////////////////////////////////////////////////////////////////////////////
//
//  Bresengham line drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//    Determine  dy=(y2-y1), dx=(x2-x1)
//    We plot first point of the line and increase errors of plotting
//    (xerr and yerr) by dx and dy respectively. If the error is greater
//    than largest from dx and dy then we should correct next point and
//    shift it to 1 pixel by that axe where error was too big.
//
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
       printf("Graphics error: %s\n", grapherrormsg(errorcode));
       printf("Press any key to halt:");
       getch();
       exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
void BLine(int x1,int y1,int x2, int y2, int color)
{
  int delta_x,delta_y,incx,incy,t,distance;
  int xerr=0,yerr=0;
  // determine dx and dy
  delta_x=x2-x1;
  delta_y=y2-y1;
  // determine steps by x and y axes (it will be +1 if we move in forward
  // direction and -1 if we move in backward direction
  if (delta_x>0) incx=1;
  else if (delta_x==0) incx=0;
  else incx=-1;
  if (delta_y>0) incy=1;
  else if (delta_y==0) incy=0;
  else incy=-1;
  delta_x=abs(delta_x);
  delta_y=abs(delta_y);
  // select largest from deltas and use it as a main distance
  if (delta_x>delta_y) distance=delta_x;
  else distance=delta_y;
  for (t=0;t<=distance+1;t++)
  {
    PutPixel(x1,y1,color);
    // increasing error
    xerr+=delta_x;
    yerr+=delta_y;
    // if error is too big then we should decrease it by changing
    // coordinates of the next plotting point to make it closer
    // to the true line
    if(xerr>distance)
    {
      xerr-=distance;
      x1+=incx;
    }
    if (yerr>distance)
    {
      yerr-=distance;
      y1+=incy;
    }
  }
}
int main(void)
{
  // initializing graphics mode
  init_gr();
  /* examples */
  BLine(0, 0, 640, 480,15);
  BLine(320, 0, 320, 480,10);
  BLine(0, 240, 640, 240,11);
  BLine(0, 480, 640, 0,12);
  BLine(320, 11, 10, 5,13);
  BLine(320, 11, 630, 5,13);
  /* clean up */
  getch();
  end_gr();
  
  return 0;
}

Рисование окружности (по Брезенхэму)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
//////////////////////////////////////////////////////////////////////////////
//
//  Bresengham circle drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//    In this algorithm all floating point math changed to sequences
//  of additions and substractions. The main idea of this algorithm
//  is to increase x and y by the value of the error between them.
//
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
float const ratio=1.0; // you can change this to draw ellipses
// This function plots points that belongs to the circle
// It recieves offsets from center for the fist quadrant
// and plots symmetrical points in all four quadrants
void plot_circle(int x,int y, int x_center, int y_center, int color)
{
  int x_start,y_start,x_end,y_end,x1,y1;
  // the distanse between start and end can be greater than 1 if ratio!=1
  y_start=y*ratio;
  y_end=(y+1)*ratio;
  x_start=x*ratio;
  x_end=(x+1)*ratio;
  for (x1=x_start;x1<x_end;++x1)
  {
    // plot points in all four quadrants
    PutPixel(x1+x_center,y+y_center,color);
    PutPixel(x1+x_center,y_center-y,color);
    PutPixel(x_center-x1,y+y_center,color);
    PutPixel(x_center-x1,y_center-y,color);
  }
  for (y1=y_start;y1<y_end;++y1)
  {
    // plot points in all four quadrants
    PutPixel(y1+x_center,x+y_center,color);
    PutPixel(y1+x_center,y_center-x,color);
    PutPixel(x_center-y1,x+y_center,color);
    PutPixel(x_center-y1,y_center-x,color);
  }
}
// This is main function that draws circle using function
void Circle(int x1,int y1,int radius, int color)
{
  int x,y,delta;
//     Y    *              we start from * and increase x step by step
//          |              decreasing y when needed
//          |
//          |
// --------------------
//          |         X
//          |
//          |
//          |
  y=radius;
  delta=3-2*radius; // delta is an error
  // calculate values for first quadrant
  for (x=0;x<y;x++) // x is a main axe
  {
    // plot points symmetrically in all quadrants
    plot_circle(x,y,x1,y1,color);
    if (delta<0) delta+=4*x+6;
    else
    {
      delta+=4*(x-y)+10;
      y--; // it's time to decrease y
    }
  }
  x=y;
  if (y!=0) plot_circle(x,y,x1,y1,color);
}
int main(void)
{
  // initializing graphics mode
  init_gr();
  /* examples */
  Circle(200,200,100,14);
  Circle(300,200,100,15);
  Circle(400,200,100,13);
  Circle(250,100,100,12);
  Circle(350,100,100,11);
  Circle(50,400,25,2);
  Circle(500,400,25,2);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Сглаживание кривой В-сплайном
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
//////////////////////////////////////////////////////////////////////////////
//
//  Curve fitting using B-splines
//  (c) Johna Smith, 1996
//
//  Method description:
//    We are drawing lines between points using the following formulas:
//    x(t)=((a3*t+a2)*t+a1)*t+a0
//    y(t)=((b3*t+b2)*t+b1)*t+b0
//    t=0..1
//    Look program for formulas for coefficients ai,bi
//    These coefficients depends on coordinates of current point,
//    previous one, next and next-next ones.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function moves CP to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void MoveTo(int x, int y)
{
  moveto(x,y);
}
// this function draws a line to (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void LineTo(int x, int y)
{
  lineto(x,y);
}
// this function draws a line from (x1,y1) to (x2,y2)
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void Line(int x1, int y1, int x2, int y2)
{
  line(x1,y1,x2,y2);
}
const N=21; // number of points
const M=300.0; // number of steps between two points
// coordinates of all given points
int x[N]={50,25,50,125,200,225,275,475,590,615,615,615,600,475,275,225,200,125,50,25,50};
int y[N]={140,100,60,40,70,90,85,75,80,85,100,115,120,125,115,110,130,160,140,100,60};
// coefficients
float t,xA,xB,xC,xD,yA,yB,yC,yD,a0,a1,a2,a3,b0,b1,b2,b3;
int n,i,j,first;
int main(void)
{
  // initializing graphics mode
  init_gr();
   /* mark the given points */
   for (i=0;i<N;i++)
   {
     Line(x[i]-4,y[i]-4,x[i]+4,y[i]+4);
     Line(x[i]+4,y[i]-4,x[i]-4,y[i]+4);
   }
   /* main loop */
   first=1;
   for(i=1;i<N-2;i++)
   {
     // calculating coefficients
     xA=x[i-1]; xB=x[i]; xC=x[i+1]; xD=x[i+2];
     yA=y[i-1]; yB=y[i]; yC=y[i+1]; yD=y[i+2];
     a3=(-xA+3*(xB-xC)+xD)/6.0;  b3=(-yA+3*(yB-yC)+yD)/6.0;
     a2=(xA-2*xB+xC)/2.0;        b2=(yA-2*yB+yC)/2.0;
     a1=(xC-xA)/2.0;             b1=(yC-yA)/2.0;
     a0=(xA+4*xB+xC)/6.0;        b0=(yA+4*yB+yC)/6.0;
     for (j=0;j<M;j++) // drawing curve between two given points
     {
       t=(float)j/(float)M;
       if (first)
       {
         first=0;
         MoveTo(((a3*t+a2)*t+a1)*t+a0,((b3*t+b2)*t+b1)*t+b0);
       } else LineTo(((a3*t+a2)*t+a1)*t+a0,((b3*t+b2)*t+b1)*t+b0);
     }
   }
   /* clean up */
   getch();
   end_gr();
   
   return 0;
}

Постепенное затемнение
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
//////////////////////////////////////////////////////////////////////////////
//
//  Fade to black (slowly turning the light off)
//  (c) Johna Smith, 1996
//
//  Method description:
//    We need only to decrease Red, Green and Blue components for each color
//    in palette. For example, if Red component of the first color was
//    initially 60 and we want to turn the light off in 30 steps we should
//    decrease Red component by 60/30=2 by step.
//  To increase productivity of this program you shouldn't use standart
//  BGI functions. Use direct methods instead.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
#include <dos.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function set filling style
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void SetFillStyle(int pattern,int color)
{
  setfillstyle(pattern,color);
}
// this function draws a bar
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void Bar(int x1, int y1, int x2, int y2)
{
  bar(x1,y1,x2,y2);
}
// this function sets palette entry using RGB color coding
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void SetRGBPalette(int colornum, int red, int green, int blue)
{
  setrgbpalette(colornum,red,green,blue);
}
// this function reads palette
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void GetPalette(struct palettetype far *palette)
{
  getpalette(palette);
}
#define N   30 // steps of fading to black
struct palettetype pal;
struct {int R;int G;int B;} colors[16];
int main(void)
{
  // initializing graphics mode
  init_gr();
  // creating palette
  GetPalette(&pal);
  for(int i=0;i<pal.size-1;i++)
  {
     colors[i+1].R=0;
     colors[i+1].G=30-2*i;
     colors[i+1].B=15+i*2;
     SetRGBPalette(pal.colors[i+1],colors[i+1].R,colors[i+1].G,colors[i+1].B);
  }
  // drawing picture
  for (i=1;i<16;i++)
  {
     SetFillStyle(SOLID_FILL,i);
     Bar((i-1)*43,0,i*43,479);
  }
  // fading to balck
  for (int j=0;j<N;j++)
  {
     //delay(50);
     for (i=1;i<pal.size;i++)
        SetRGBPalette(pal.colors[i],
   colors[i].R*(1-(float)j/N),
          colors[i].G*(1-(float)j/N),
          colors[i].B*(1-(float)j/N));
  }
  // and back to the light
  for (j=0;j<N;j++)
  {
     //delay(50);
     for (i=1;i<pal.size;i++)
        SetRGBPalette(pal.colors[i],
          colors[i].R*(float)j/N,
          colors[i].G*(float)j/N,
          colors[i].B*(float)j/N);
  }
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Заливка замкнутой области
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
//////////////////////////////////////////////////////////////////////////////
//
//  Filling closed area using given color
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Fill all pixels on the current line
//    2) Call FloodFill recursively for two remote points to reduce
//       the depth of the recursion
//    3) Call FloodFill recursively for all points on two nearest
//       horizontal lines
//  ! This program is NOT optimized for best performance !
//  To do so don't use PutPixel & GetPixel function -
//  change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
// this function gets color of pixel in (x,y) position
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
int GetPixel(int x, int y)
{
  return getpixel(x,y);
}
// this function gets maximum horizontal coordinate
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
int GetMaxX(void)
{
  return getmaxx();
}
// this function gets maximum vertical coordinate
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
int GetMaxY(void)
{
  return getmaxy();
}
void FloodFill(int x, int y,int fill_color, int border_color)
{
  int x_left,x_right,YY,i;
  int XMAX,YMAX;
  XMAX=GetMaxX();
  YMAX=GetMaxY();
  // Light as many pixels as possible on line y
  // and determine x_left and x_right
  x_left=x_right=x;
  while (GetPixel(x_left,y)!=border_color && x_left>0)
  {
     PutPixel(x_left,y,fill_color);
     x_left--;
  }
  x_right++;
  while (GetPixel(x_right,y)!=border_color && x_right<XMAX)
  {
     PutPixel(x_right,y,fill_color);
     x_right++;
  }
  // Recursive calls of FloodFill for two remote points
  x=(x_right+x_left)>>1; //shifting means division by 2
  for(i=-1;i<=1;i+=2)
  {
     YY=y;
     while (GetPixel(x,YY)!=border_color && YY<YMAX && YY>0) YY+=i;
     YY=(y+YY)>>1;
     if (GetPixel(x,YY)!=border_color && GetPixel(x,YY)!=fill_color) FloodFill(x,YY,fill_color,border_color);
  }
  // Recursive calls for all "dark" (not filled) pixels next
  // to the line y (with x values between x_left and x_right
  for(YY=y-1;YY<=y+1;YY+=2)
  {
    x=x_left+1;
    while (x<x_right && YY>0 && YY<YMAX)
    {
      if (GetPixel(x,YY)!=border_color && GetPixel(x,YY)!=fill_color) FloodFill(x,YY,fill_color,border_color);
      x++;
    }
  }
}
int main(void)
{
   // initializing graphics mode
   init_gr();
   // drawing area to fill
   moveto(520,90); lineto(520,190); lineto(470,240);
   lineto(520,290); lineto(520,390);
   lineto(320,315); lineto(120,390);
   lineto(120,290); lineto(170,240);
   lineto(120,190); lineto(120, 90);
   lineto(270,90); lineto(270,120);
   lineto(170,120); lineto(170,150);
   lineto(470,150); lineto(470,120);
   lineto(370,120); lineto(370,90); lineto(520,90);
   moveto(370,240); lineto(320,290); lineto(270,240);
   lineto(320,190); lineto(370,240);
   /* example */
   FloodFill(121,91,1,15);
   FloodFill(320,240,9,15);
   FloodFill(100,20,3,15);
   /* clean up */
   getch();
   end_gr();
   return 0;
}

Рисование линии
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
//////////////////////////////////////////////////////////////////////////////
//
//  Line drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//    Determine k=dy/dx, dy=(y2-y1), dx=(x2-x1)
//    Line equation is: y(x)=y1+x*k so we need only to plot all
//    these points. This algorithm is very slow, because it use
//    floating point math (requires coprocessor for best productivity)
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
void Line(int x1,int y1,int x2, int y2, int color)
{
  int delta_x,delta_y,incx,incy;
  float k;
  // determine dx and dy
  delta_x=x2-x1;
  delta_y=y2-y1;
  // determine steps by x and y axes (it will be +1 if we move in forward
  // direction and -1 if we move in backward direction
  if (delta_x>0) incx=1;
  else if (delta_x==0) incx=0;
  else incx=-1;
  if (delta_y>0) incy=1;
  else if (delta_y==0) incy=0;
  else incy=-1;
  delta_x=abs(delta_x);
  delta_y=abs(delta_y);
  // select greatest from deltas and use it as a main axe
  if (delta_x>delta_y && delta_x!=0)
  {
     k=(float)delta_y/delta_x;
     for (int i=0;i<delta_x;i++) PutPixel(x1+incx*i,y1+incy*floor(k*i),color);
  }
  else
  {
     k=(float)delta_x/delta_y;
     for (int i=0;i<delta_y;i++) PutPixel(x1+incx*floor(k*i),y1+incy*i,color);
  }
}
int main(void)
{
  // initializing graphics mode
  init_gr();
  /* examples */
  Line(0, 0, 640, 480,15);
  Line(320, 0, 320, 480,10);
  Line(0, 240, 640, 240,11);
  Line(0, 480, 640, 0,12);
  Line(320, 11, 10, 5,13);
  Line(320, 11, 630, 5,13);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Рисование прямоугольника
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
//////////////////////////////////////////////////////////////////////////////
//
//  Rectangle drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//  Draw four lines, using simplified line drawing method, cause we know
//  that only vertical and horizontal lines will be drawn.
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
void Line(int x1,int y1,int x2, int y2, int color)
{
  int delta_x,delta_y,incx,incy;
  // determine dx and dy
  delta_x=x2-x1;
  delta_y=y2-y1;
  // determine steps by x and y axes (it will be +1 if we move in forward
  // direction and -1 if we move in backward direction
  if (delta_x>0) incx=1;
  else if (delta_x==0) incx=0;
  else incx=-1;
  if (delta_y>0) incy=1;
  else if (delta_y==0) incy=0;
  else incy=-1;
  delta_x=abs(delta_x);
  delta_y=abs(delta_y);
  // select greatest from deltas and use it as a main axe
  if (delta_x!=0)
     for (int i=0;i<delta_x;i++) PutPixel(x1+incx*i,y1,color);
  else
     for (int i=0;i<delta_y;i++) PutPixel(x1,y1+incy*i,color);
}
// this function draws a rectangle
void Rectangle(int x1, int y1, int x2, int y2, int color)
{
  Line(x1,y1,x1,y2,color);
  Line(x1,y2,x2,y2,color);
  Line(x2,y2,x2,y1,color);
  Line(x2,y1,x1,y1,color);
}
int main(void)
{
  // initializing graphics mode
  init_gr();
  /* examples */
  Rectangle(0,0,639,479,10);
  Rectangle(100,200,300,400,11);
  Rectangle(600,100,20,11,12);
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Фрактал (листья папоротника)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
//////////////////////////////////////////////////////////////////////////////
//
//  Leafs drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//    We take a closed area, which will be main leaf. We want to generate
//  four fractal leafs from this one. To make it we need to scale, move and
//  rotate main leaf. These transformations described by the following
//  equations:
//          x'=ax+by+c, y'=dx+ey+f
//    So we need 6 coefficients for each leaf (24 coefficients at all).
//  These coefficients are calculated with special mathematical method
//  (here are already calculated coefficients)
//    To draw whole picture we need an iterational process:
//   1) Take any point within main leaf borders
//   2) Call function iterate for one of 4 leafs (its nuber we select
//      randomly but than the larger is leaf's area the higher is probability
//      that we call 'iterate' for this leaf
//   3) Then we take transformed point as (x,y) and repeat iterations
//
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <stdlib.h>
#include <dos.h>
#include <graphics.h>
#include <conio.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
    /* request autodetection */
    int gdriver = DETECT, gmode, errorcode;
    /* initialize graphics mode */
    initgraph(&gdriver, &gmode, "");
    /* read result of initialization */
    errorcode = graphresult();
    if (errorcode != grOk)    /* an error occurred */
    {
       printf("Graphics error: %s\n", grapherrormsg(errorcode));
       printf("Press any key to halt:");
       getch();
       exit(1);               /* return with error code */
    }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
#define N       200000L // number of points
float x=0, y=0; // iteration variables (are global)
//                       a     b     c    d    e    f        probability
float coeff[4][6] = {{ 0.00, 0.00, 0.00,0.16,0.00,0.00,},   // 01%
                     { 0.85, 0.04,-0.04,0.85,0.00,1.60,},   // 85%
                     { 0.05,-0.26, 0.23,0.22,0.00,1.60,},   // 07%
                     {-0.15, 0.30, 0.26,0.24,0.00,0.44,},}; // 07%
int colors[5]={10,11,2,3,2};
void iterate(char i,int c)
{
  float x1,y1;
  x1=x*coeff[i][0]+y*coeff[i][1]+coeff[i][4];
  y1=x*coeff[i][2]+y*coeff[i][3]+coeff[i][5];
  x=x1;
  y=y1;
  putpixel ((int)(y*64),240-(int)(x*48),c);
}
int main (void)
{
  int leaf,color;
  // initializing graphics mode
  init_gr();
  // drawing leafs
  for (long int i=0;i<N;i++)
  {
     color=colors[random(5)];  // random color
     leaf=random(1000);  // selecting leaf to draw
     if (leaf<=10) iterate(0,color);
     else if (leaf<=860) iterate(1,color);
     else if (leaf<=930) iterate(2,color);
     else iterate(3,color);
  }
  /* clean up */
  getch();
  end_gr();
  return 0;
}

Множество Мандельброта
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
//////////////////////////////////////////////////////////////////////////////
//
//  Mandelbrot set drawing
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Assume z(0)=0
//    2) Organize two loops by X and Y axes
//    3) Determine constant C=X+iY
//    4) Using iterative formula z(i)=z(i-1)*z(i-1)+c calculate Z(DEPTH),
//       where DEPTH is max depth of iterations
//    5) If Z(DEPTH) is near Z(0) then point (X,Y) belongs to Mandelbrot set
//
//  ! This program is NOT optimized for best performance !
//  To do so don't use putpixel function - change bytes in video memory directly.
//
//////////////////////////////////////////////////////////////////////////////
#define DEPTH 100 // iterations depth
#include <graphics.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include <math.h>
// this function initializes graphics mode
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void init_gr(void)
{
   /* request autodetection */
   int gdriver = DETECT, gmode, errorcode;
   /* initialize graphics mode */
   initgraph(&gdriver, &gmode, "");
   /* read result of initialization */
   errorcode = graphresult();
   if (errorcode != grOk)    /* an error occurred */
   {
      printf("Graphics error: %s\n", grapherrormsg(errorcode));
      printf("Press any key to halt:");
      getch();
      exit(1);               /* return with error code */
   }
}
// this function shuts graphics mode down
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void end_gr(void)
{
  closegraph();
}
// this function puts pixel on the screen in (x,y) position using color 'color'
// it will work only if you're using Borland C++ compiler & BGI drivers
// if you're using another compiler you should overwrite body of this function
void PutPixel(int x, int y, int color)
{
  putpixel(x,y,color);
}
void Mandelbrot(void)
{
  float zi, zr, ci, cr, tmp;
  // assuming graphics mode 640x480 16 colors
  for(int i=-320;i<320;i++)  // for all pixels on the X axe
  {
     ci=((float)i)/320.0; // setting Im c = i/320
     for(int j=-380;j<160;j++)  // for all pixels on the Y axe
     {
        cr=((float)j)/240.0; // setting Re c = j/240
        zi=zr=0.0; // fixing z=0 we'll change c - it's Mandelbrot set
        for(int k=0;k<DEPTH;k++)
        {
   // z=z*z+c
   //(zr+i*zi)*(zr+i*zi)=(zr*zr-zi*zi)+i*(2*zr*zi)
   // zi=zr*zr-zi*zi+cr
   //
   tmp=zr*zr-zi*zi;
   zi=2*zr*zi+ci;
   zr=tmp+cr;
   if (zr*zr+zi*zi>1.0E16) break; // break if |z| is very big
    }
    if (k<DEPTH)
   PutPixel(i+320,j+380,k%3+1); // z was very big => it is external point
    else PutPixel(i+320,j+380,11); // internal point
     }
     if(kbhit()) break;
  }
}
int main(void)
{
  // initializing graphics mode
  init_gr();
  /* example */
  Mandelbrot();
  /* clean up */
  getch();
  end_gr();
  return 0;
}
Лучшие ответы (1)
После регистрации реклама в сообщениях будет скрыта и будут доступны все возможности форума.
LK
Заблокирован
19.05.2013, 15:15  [ТС]     Коллекция алгоритмов от Johna Smith #2
Специальные функции

arccos x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  arccos x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arccos x using the following formula:
//                                        2n+1
//              Pi       N  1*3*5..(2n-1)x
//   arccos x = - - x - SUM ------------------ , |x|<1
//              2       n=1 2*4*6..(2n)(2n+1)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi      3.1415926536
#define N       30
double Arccos(double x)
{
  double arccos=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x*(2*i-1)*(2*i-1)/((2*i+1)*2*i);
    arccos+=a;
  }
  return Pi/2-arccos;
}
void main(void)
{
  printf("C++ arccos(0.5)=%g, this arccos(0.5)=%g.",acos(0.5),Arccos(0.5));
}

arccos z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function arccos z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//       1            2    2     1            2    2
//    A= - sqrt( (x+1)  + y  ) + - sqrt( (x-1)  + y  )
//       2                       2
//
//       1            2    2     1            2    2
//    B= - sqrt( (x+1)  + y  ) - - sqrt( (x-1)  + y  )
//       2                       2
//                                       2
//    arccos z = acrcos B - i*ln(a+sqrt(a  -1))
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arccos(complex z)
{
  complex c;
  float A,B;
  A=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 +
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  B=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 -
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  c.re=acos(B);
  c.im=-log(A+sqrt(A*A-1));
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Arccos(");
  show_complex(z);
  printf(") = ");
  show_complex(Arccos(z));
}

arch x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  arcch x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arccos x using the following formula:
//
//                       N  1*3*5..(2n-1)       1
//   arcch x = ln(2x) - SUM --------------- * ------  , x>1
//                      n=1 2*4*6..(2n)(2n)   x^(2n)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Arch(double x)
{
  double arch=1/(2*x*2*x);
  double a=1/(2*x*2*x);
  for (int i=2;i<N;i++)
  {
    a *= (2*i-1)*(2*i-2)/(2*i*2*i*x*x);
    arch += a;
  }
  return log(2*x)-arch;
}
void main(void)
{
  printf("C++ arch(6.13229)=2.5, this arccos(6.13229)=%g.",Arch(6.13229));
}

arch z
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function arccos z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//       1            2    2     1            2    2
//    A= - sqrt( (x+1)  + y  ) + - sqrt( (x-1)  + y  )
//       2                       2
//
//       1            2    2     1            2    2
//    B= - sqrt( (x+1)  + y  ) - - sqrt( (x-1)  + y  )
//       2                       2
//                            2
//    arch z = +/- ln(a+sqrt(a  -1)) +/- i*acrcos B
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arch(complex z)
{
  complex c;
  float A,B;
  A=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 +
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  B=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 -
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  c.re=log(A+sqrt(A*A-1));
  c.im=-acos(B);
  return c;
}
complex z={3,2};
void main(void)
{
  complex c;
  printf("Arch(");
  show_complex(z);
  c=Arch(z);
  printf(") = +/- %f +/-i%f",c.re,c.im);
}

arcsin x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  arcsin x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arcsin x using the following formula:
//                                      2n+1
//                     N  1*3*5..(2n-1)x
//    arcsin x  = x + SUM ------------------, |x|<1
//                    n=1 2*4*6..(2n)(2n+1)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Arcsin(double x)
{
  double arcsin=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a *= x*x*(2*i-1)*(2*i-1)/((2*i+1)*2*i);
    arcsin += a;
  }
  return arcsin;
}
void main(void)
{
  printf("C++ arcsin(0.5)=%g, this arcsin(0.5)=%g.",asin(0.5),Arcsin(0.5));
}

arcsin z
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function arcsin z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//       1            2    2     1            2    2
//    A= - sqrt( (x+1)  + y  ) + - sqrt( (x-1)  + y  )
//       2                       2
//
//       1            2    2     1            2    2
//    B= - sqrt( (x+1)  + y  ) - - sqrt( (x-1)  + y  )
//       2                       2
//                                       2
//    arcsin z = acrsin B + i*ln(a+sqrt(a  -1))
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arcsin(complex z)
{
  complex c;
  float A,B;
  A=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 +
  sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  B=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 -
  sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  c.re=asin(B);
  c.im=log(A+sqrt(A*A-1));
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Arcsin(");
  show_complex(z);
  printf(") = ");
  show_complex(Arcsin(z));
}

arctg x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
//////////////////////////////////////////////////////////////////////////////
//
//  arctg x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arctg x using the following formulas:
//                              2n+1
//                  N      n   x
//     arctg x  =  SUM (-1)  --------, |x|<1
//                 n=0         2n+1
//
//
//                  N      n+1       1
//     arctg x  =  SUM (-1)    --------------, |x|>1
//                 n=0         (2n+1)x^(2n+1)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       100
#define Pi      3.1415926536
double Arctg(double x)
{
  double arctg;
  double a;
  if (fabs(x)<1)
  {
    a=x;
    arctg=x;
    for (int i=1;i<N;i++)
    {
       a *= -x*x;
       arctg += a/(2*i+1);
    }
  } 
  else
  {
    a=-1/x;
    arctg=-1/x;
    for (int i=1;i<N;i++)
    {
       a *= -1/(x*x);
       arctg += a/(2*i+1);
    }
    arctg += (x>0?Pi/2.0:-Pi/2.0);
  }
  return arctg;
}
void main(void)
{
  printf("C++ arctg(0.5)=%g, this arctg(0.5)=%g.",atan(0.5),Arctg(0.5));
}

arctg z
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
//////////////////////////////////////////////////////////////////////////////
//
//  arctg z function calculating, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//                 1          2x          1    x^2+(y+1)^2
//     arctg z  =  - arctg(---------) + i - ln(-----------)
//                 2       1-x^2-y^2      4    x^2+(y-1)^2
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arctg(complex z)
{
  complex c;
  
  c.re=atan(2*z.re/(1-z.re*z.re-z.im*z.im))/2;
  c.im=log((z.re*z.re+(z.im+1)*(z.im+1))/(z.re*z.re+(z.im-1)*(z.im-1)))/4;
  return c;
}
complex z={0.3,0.2};
void main(void)
{
  printf("Arctg(");
  show_complex(z);
  printf(") = ");
  show_complex(Arctg(z));
}

arcsh x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  arsh x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arccos x using the following formula:
//                                          2n+1
//                   N      n 1*3*5..(2n-1)x
//     arsh x = x + SUM (-1)  ------------------ , |x|<1
//                  n=1       2*4*6..(2n)(2n+1)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Arsh(double x)
{
  double arsh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a *= -x*x*(2*i-1)*(2*i-1)/((2*i+1)*2*i);
    arsh += a;
  }
  return arsh;
}
void main(void)
{
  printf("Exact arsh(0.5210953)=0.5, this arsh(0.5)=%g.",Arsh(0.5210953));
}

arcsh z
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function arcsh z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//       1            2    2     1            2    2
//    A= - sqrt( (x+1)  + y  ) + - sqrt( (x-1)  + y  )
//       2                       2
//
//       1            2    2     1            2    2
//    B= - sqrt( (x+1)  + y  ) - - sqrt( (x-1)  + y  )
//       2                       2
//                         2
//    arcsh z = ln(a+sqrt(a  -1)) - i*acrsin B
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arcsh(complex z)
{
  complex c;
  float A,B,swp;
  swp=z.im;    // calculating arsh(z) as -arcsin(iz)
  z.im=z.re;
  z.re=-swp;
  A=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 +
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  B=sqrt((z.re+1)*(z.re+1)+z.im*z.im)/2 -
         sqrt((z.re-1)*(z.re-1)+z.im*z.im)/2;
  c.re=log(A+sqrt(A*A-1));
  c.im=-asin(B);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Arcsh(");
  show_complex(z);
  printf(") = ");
  show_complex(Arcsh(z));
}

arth x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  arth x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating arth x using the following formulas:
//                       2n+1
//                  N   x
//      arth x  =  SUM -------
//                 n=0  2n+1
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Arth(double x)
{
  double arth=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a *= x*x;
    arth += a/(2*i+1);
  }
  return arth;
}
void main(void)
{
  printf("Exact arth(0.4621171)=0.5, this arth(0.4621171)=%g.",Arth(0.4621171));
}

arth z
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
//////////////////////////////////////////////////////////////////////////////
//
//  arctg z function calculating, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//                   1    x^2+(y+1)^2    1          2x
//     arctg z  =    - ln(-----------) -i- arctg(---------)
//                   4    x^2+(y-1)^2    2       1-x^2-y^2
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Arth(complex z)
{
  complex c;
  float swp;
  swp=z.im;    // calculating arth(z) as -arctg(iz)
  z.im=z.re;
  z.re=-swp;
  c.re=log((z.re*z.re+(z.im+1)*(z.im+1))/(z.re*z.re+(z.im-1)*(z.im-1)))/4;
  c.im=-atan(2*z.re/(1-z.re*z.re-z.im*z.im))/2;
  return c;
}
complex z={0.3,0.2};
void main(void)
{
  printf("Arth(");
  show_complex(z);
  printf(") = ");
  show_complex(Arth(z));
}

Функция Бесселя
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
//////////////////////////////////////////////////////////////////////////////
//
//  Bessel functions calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Bessel functions are solutions of the followind eqation:
//         2                     2
//        d w     1 dw          v
//        ---- +  - -- + ( 1 - --- )w = 0
//          2     x dx           2
//        dx                    x
//
//    This program calculates Bessel functions for integer v as infinite sum
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
// this function returns value of Jn(x)
double Jn(int n,double x)
{
  float nfact=1,ifact=1,nifact=1,sum=0,x2,element;
  int i;
  if (n>1) for (i=0;i<n;i++) nfact*=(i+1);  // n factorial
  x2=x*x/4;
  i=0;
  do
  {
    i++;
    ifact*=i;
    nifact*=(n+i);
    element=(i%2==0?1:-1)*pow(x2,i)/(ifact*nifact);
    sum+=element;
  } while (fabs(element)>1e-9);
  return pow(x/2,n)*(sum+1)/nfact;
}
// this function returns value of In(x)
double In(int n, double x)
{
  float nfact=1,ifact=1,nifact=1,sum=0,x2,element;
  int i;
  if (n>1) for (i=0;i<n;i++) nfact*=(i+1);  // n factorial
  x2=x*x/4;
  i=0;
  do
  {
    i++;
    ifact*=i;
    nifact*=(n+i);
    element=pow(x2,i)/(ifact*nifact);
    sum+=element;
  } while (fabs(element)>1e-9);
  return pow(x/2,n)*(sum+1)/nfact;
}
void main(void)
{
  // Examples
  printf("J0(0.5)=%f\n",Jn(0,0.5));
  printf("J30(20)=%f\n",Jn(30,20));
  printf("I0(2)=%f\n",In(0,2));
  printf("I1(2)=%f\n",In(1,2));
}

Гамма-функция
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
//////////////////////////////////////////////////////////////////////////////
//
//  Gamma function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//             infinity  -t  x-1
//   Gamma(x)= Integral e   t    dt
//                0
//
//   This integral is approximated and calculated by Stirling's formula:
//                A
//             (A)    1/12A
//      Gamma= (-)  e       D sqrt(2 Pi A)
//             (e)
//
//   This algorithm is oriented for -10<x<=70
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define E  2.7182818285
#define Pi  3.1415926536
// this function returns value of Gamma(x)
double Gamma(double x)
{
  double A,z,g,c,D;
  int b;
  A=x-1;
  z=A-64;
  if (z<0)
  {
    b=abs((int)(z-1));
    c=A;
    A+=b;
    D=1;
    for(int i=b;i>0;i--)
    {
       c++;
       D/=c;
    }
  } else D=1;
  g=pow(A/E,A)*exp(1/(12*A))*D*sqrt(2*Pi*A);
  return g;
}
void main(void)
{
  // Examples
  printf("Gamma(-1.5)=%f, exact value is 2.363271801\n",Gamma(-1.5));
  printf("Gamma(-0.5)=%f, exact value is -3.544907702\n",Gamma(-0.5));
  printf("Gamma(1.5)=%f, exact value is 0.886226926\n",Gamma(1.5));
  printf("Gamma(40)=%f, exact value is 2.039788208E+046\n",Gamma(40));
}

Бета-функция
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
//////////////////////////////////////////////////////////////////////////////
//
//  Beta function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//                 1      x-1      w-1
//   Beta(x,w)= Integral t    (1-t)   dt
//                 0
//              Gamma(x) Gamma(w)
//   Beta(x,w)= -----------------
//                 Gamma(x+w)
//
//   This algorithm is oriented for -10 < x,w,x+w <= 70
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define E       2.7182818285
#define Pi      3.1415926536
// this function returns value of Gamma(x)
double Gamma(double x)
{
  double A,z,g,c,D;
  int b;
  A=x-1;
  z=A-64;
  if (z<0)
  {
     b=abs((int)(z-1));
     c=A;
     A+=b;
     D=1;
     for(int i=b;i>0;i--)
     {
        c++;
        D/=c;
     }
  } else D=1;
  g=pow(A/E,A)*exp(1/(12*A))*D*sqrt(2*Pi*A);
  return g;
}
double Beta(double x, double w)
{
  return(Gamma(x)*Gamma(w)/Gamma(x+w));
}
void main(void)
{
  // Examples
  printf("Beta(1,2)=%f, exact value is 0.5\n",Beta(1,2));
}

ch x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  ch x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating sh x using the following formulas:
//                      2n
//                  N  x
//        ch x  =  SUM ---
//                 n=0 2n!
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a *= x*x/((2*i+2)*(2*i+1));
    ch += a;
  }
  return ch;
}
void main(void)
{
  printf("Exact ch(0.5)=%g, this ch(0.5)=%g.",(exp(0.5)+exp(-0.5))/2,Ch(0.5));
}

ch z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function Ch(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  ch z = ch x cos y + ish x sin y
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex complexCh(complex z)
{
  complex c;
  c.re=Ch(z.re)*cos(z.im);
  c.im=Sh(z.re)*sin(z.im);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Ch(");
  show_complex(z);
  printf(") = ");
  show_complex(complexCh(z));
}

cosec x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Cosecans function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating cosec x using the following formula:
//
//            1   1     7    3    31    5    127    7
//  cosec x = - + - x + --- x  + ----- x  + ------ x  + ...
//            x   6     360      15120      604800
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Cosec(double x)
{
  double cosec=0;
  double coefficient[]={1, 1/6.0, 7/360.0, 31/15120.0, 127/604800.0};
  double a=1/x;
  for (int i=0;i<sizeof(coefficient)/sizeof(double);i++)
  {
    cosec+=a*coefficient[i];
    a*=x*x;
  }
  return cosec;
}
void main(void)
{
  printf("C++ cosec(0.5)=%g, this cosec(0.5)=%g.",1/sin(0.5),Cosec(0.5));
}

cos x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Cosine function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating cos x using the following formula:
//                  2     4     6     8
//                 x     x     x     x
//    cos x = 1 - --- + --- - --- + --- - ...
//                 2!    4!    6!    8!
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       300
double Cos(double x)
{
  double cos=0;
  double a=1;
  for (int i=1;i<N;i++)
  {
    cos+=a;
    a=-a*x*x/((2*i-1)*2*i);
  }
  return cos;
}
void main(void)
{
  printf("C++ Cos(0.5)=%g, this Cos(0.5)=%g.",cos(0.5),Cos(0.5));
}

cos z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function cos(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  cos z = cos x ch y + isin x sh y
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Cos(complex z)
{
  complex c;
  c.re=cos(z.re)*Ch(z.im);
  c.im=sin(z.re)*Sh(z.im);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Cos(");
  show_complex(z);
  printf(") = ");
  show_complex(Cos(z));
}

ctg x
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Cotangent function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating ctg x using the following formula:
//
//            1   1     1   3    2   5    1    7
//    ctg x = - - - x - -- x  - --- x  - ---- x  - ...
//            x   3     45      945      4725
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Ctg(double x)
{
  double ctg=0;
  double coefficient[]={1, -1/3.0, -1/45.0, -2/945.0, -1/4725.0};
  double a=1/x;
  for (int i=0;i<sizeof(coefficient)/sizeof(double);i++)
  {
    ctg += a*coefficient[i];
    a *= x*x;
  }
  return ctg;
}
void main(void)
{
  printf("C++ ctg(0.5)=%g, this ctg(0.5)=%g.",1/tan(0.5),Ctg(0.5));
}

e^x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  e^x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating e^x using the following formula:
//                   2    3
//     x       x    x    x  
//    e  = 1 + -- + -- + -- + ...
//             1!   2!   3!
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       300
double Exp(double x)
{
  double exp=0;
  double a=1;
  for (int i=1;i<N;i++)
  {
    exp+=a;
    a*=x/i;
  }
  return exp;
}
void main(void)
{
  printf("C++ exp(10)=%g, this exp(10)=%g.",exp(10),exp(10));
}

e^z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function e^z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  e^z = e^x cos y + i e^x sin y
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Exp(complex z)
{
  complex c;
  c.re=exp(z.re)*cos(z.im);
  c.im=exp(z.re)*sin(z.im);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Exp(");
  show_complex(z);
  printf(") = ");
  show_complex(Exp(z));
}

tg x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Tangent function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating tg x using the following formula:
//
//                1  3   2   5   17   7   62    9
//     tg x = x + - x  + -- x  + --- x  + ---- x  + ...
//                3      15      315      2835
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Tg(double x)
{
  double tg=0;
  double coefficient[]={1,1/3.0,2/15.0,17/315.0,62/2835.0};
  double a=x;
  for (int i=0;i<sizeof(coefficient)/sizeof(double);i++)
  {
    tg+=a*coefficient[i];
    a*=x*x;
  }
  return tg;
}
void main(void)
{
  printf("C++ Tg(0.5)=%g, this Tg(0.5)=%g.",tan(0.5),Tg(0.5));
}

tg z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function tg(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//             sin(2x)             sh(2y)
//  tg z = --------------- + i --------------
//         cos(2x)+ch(2*y)     cos(2x)+ch(2y)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Tg(complex z)
{
  complex c;
  c.re=sin(2*z.re)/(cos(2*z.re)+Ch(2*z.im));
  c.im=Sh(2*z.im)/(cos(2*z.re)+Ch(2*z.im));
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Tg(");
  show_complex(z);
  printf(") = ");
  show_complex(Tg(z));
}

th x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  th x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating th x using the following formula:
//
//                1  3   2   5   17   7    62   9                Ї
//    sec x = x - - x  + -- x  - --- x  + ---- x  - ... ,  |x| < -
//                3      15      315      2835                   2
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Th(double x)
{
  double th=0;
  double coefficient[]={1, -1/3.0, 2/15.0, -17/315.0, 62/2835.0};
  double a=x;
  for (int i=0;i<sizeof(coefficient)/sizeof(double);i++)
  {
    th += a*coefficient[i];
    a *= x*x;
  }
  return th;
}
void main(void)
{
  printf("Exact th(0.5)=%g, this th(0.5)=%g.",
         (exp(0.5)-exp(-0.5))/(exp(0.5)+exp(-0.5)),Th(0.5));
}

th z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function th(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//              sh(2x)             sin(2y)
//  th z = --------------- + i --------------
//         ch(2x)+cos(2*y)     ch(2x)+cos(2y)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex complexTh(complex z)
{
  complex c;
  c.re=Sh(2*z.re)/(Ch(2*z.re)+cos(2*z.im));
  c.im=sin(2*z.im)/(Ch(2*z.re)+cos(2*z.im));
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Th(");
  show_complex(z);
  printf(") = ");
  show_complex(complexTh(z));
}

sin x (x- real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Sine function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating sin x using the following formula:
//                  3     5     7     9
//                 x     x     x     x
//    sin x = x - --- + --- - --- + --- - ...
//                 3!    5!    7!    9!
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       300
double Sin(double x)
{
  double sin=0;
  double a=x;
  for (int i=1;i<N;i++)
  {
    sin+=a;
    a=-a*x*x/((2*i+1)*2*i);
  }
  return sin;
}
void main(void)
{
  printf("C++ Sin(0.5)=%g, this Sin(0.5)=%g.",sin(0.5),Sin(0.5));
}

sin z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function sin(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  sin z = sin x ch y + icos x sh y
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Sin(complex z)
{
  complex c;
  c.re=sin(z.re)*Ch(z.im);
  c.im=cos(z.re)*Sh(z.im);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Sin(");
  show_complex(z);
  printf(") = ");
  show_complex(Sin(z));
}

sec x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Secans function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating sec x using the following formula:
//
//                1  2   5   4   61   6   277   8
//    sec x = 1 + - x  + -- x  + --- x  + ---- x  + ...
//                2      24      720      8064
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Sec(double x)
{
  double sec=0;
  double coefficient[]={1, 1/2.0, 5/24.0, 61/720.0, 277/8064.0};
  double a=1;
  for (int i=0;i<sizeof(coefficient)/sizeof(double);i++)
  {
    sec+=a*coefficient[i];
    a*=x*x;
  }
  return sec;
}
void main(void)
{
  printf("C++ Sec(0.5)=%g, this Sec(0.5)=%g.",1/cos(0.5),Sec(0.5));
}

sh x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  sh x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating sh x using the following formulas:
//                       2n-1
//                  N   x
//        sh x  =  SUM -------
//                 n=1 (2n-1)!
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       30
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void main(void)
{
  printf("Exact sh(0.5)=%g, this sh(0.5)=%g.",(exp(0.5)-exp(-0.5))/2,Sh(0.5));
}

sh z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function Sh(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  sh z = sh x cos y + ich x sin y
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
#define N       30  // precision
double Ch(double x)
{
  double ch=1;
  double a=1;
  for (int i=0;i<N;i++)
  {
    a*=x*x/((2*i+2)*(2*i+1));
    ch+=a;
  }
  return ch;
}
double Sh(double x)
{
  double sh=x;
  double a=x;
  for (int i=1;i<N;i++)
  {
    a*=x*x/((2*i)*(2*i+1));
    sh+=a;
  }
  return sh;
}
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex complexSh(complex z)
{
  complex c;
  c.re=Sh(z.re)*cos(z.im);
  c.im=Ch(z.re)*sin(z.im);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Sh(");
  show_complex(z);
  printf(") = ");
  show_complex(complexSh(z));
}

Факториал (формула Стирлинга)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Factorial calculation (Stirling solution)
//  (c) Johna Smith, 1996
//
//  Method description:
//   We just use Stirling method, improving it by using fast method
//  of calculating n^n.
//   !! Warning! this formula gives best results ONLY if n>>10
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi 3.1415926536
long double factorial (long int n)
{
  long double factorial,b=1,tmp;
  long int k;
  factorial=sqrt(2*Pi*n)*exp(-n+(1-1/(30*n))/(12*n));
  // quick powering n^n
  tmp=n;
  k=n;
  while (k!=0)
  {
    if (k%2!=0) b*=tmp;
    k/=2;
    tmp*=tmp;
  }
  factorial*=b;
  return factorial;
}
void main(void)
{
  printf("5! = %Lf\n35! = %Lf\n100! = %Lf\n",factorial(5L),
                        factorial(35L),factorial(100L));
}

sqrt(z), where z is a complex value
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function sqrt(z), where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//
//  sqrt(z) = +/- sqrt((x+sqrt(x^2+y^2))/2) +/- i*sqrt(-x+sqrt(x^2+y^2))/2)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Sqrt(complex z)
{
  complex c;
  c.re=sqrt((z.re+sqrt(z.re*z.re+z.im*z.im))/2);
  c.im=sqrt((-z.re+sqrt(z.re*z.re+z.im*z.im))/2);
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Sqrt(");
  show_complex(z);
  z=Sqrt(z);
  printf(") = +/- %f +/- i*%f",z.re,z.im);
}

ln x (x - real)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  ln x function calculating.
//  (c) Johna Smith, 1996
//
//  Method description:
//    Calculating ln x using the following formula:
//
//               N     (x-1)^(2n+1)
//    ln x  = 2 SUM ------------------
//              n=0 (2n+1)(x+1)^(2n+1)
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       300
double Ln(double x)
{
  double ln=0;
  double a=(x+1)/(x-1);
  for (int i=0;i<N;i++)
  {
    a*=(x-1)*(x-1)/((x+1)*(x+1));
    ln+=a/(2*i+1);
  }
  return 2*ln;
}
void main(void)
{
  printf("C++ ln(10)=%g, this ln(10)=%g.",log(10),Ln(10));
}

ln z (z - complex)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Calculactin function ln z, where z is a complex value
//  (c) Johna Smith, 1996
//
//  Method description:
//         1   2   2             y
//  ln z = - (x + y ) + i*(arctg - +2kЇ), k=0, +/-1, +/-2, ...
//         2                     x
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Ln(complex z)
{
  complex c;
  c.re=log(z.re*z.re+z.im*z.im)/2;
  c.im=atan(z.im/z.re); // you can add 2kЇ here (we assume that k=0)
  return c;
}
complex z={3,2};
void main(void)
{
  printf("Ln(");
  show_complex(z);
  printf(") = ");
  show_complex(Ln(z));
}


Динамические структуры данных

Bidirectional cycled list
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
//////////////////////////////////////////////////////////////////////////////
//
//  Dynamic structures (bidirectional cycled list)
//  (c) Johna Smith, 1996
//
//  Method description:
// .................
// -> * -> * -> * ->
// <- * <- * <- * <-
// .  A    B    C  .
// .................
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <alloc.h>
struct item
{
  int element;
  item *prev;
  item *next;
};
item *list; // base element of the list
// this function removes element pointed by i from list
void Remove(item* i)
{
  i->next->prev=i->prev;
  i->prev->next=i->next;
  free(i);
}
// this function inserts element after
// element pointed by i
void Insert(item *i, int element)
{
  item *tmp;
  // creating new item
  tmp=(item*)malloc(sizeof(item));
  tmp->element=element;
  tmp->next=i->next;
  tmp->prev=i->next->prev;
  // correcting nearest elements pointers
  i->next->prev=tmp;
  i->next=tmp;
}
// this function searches element in the list
item* Search(int element)
{
  item *p,*q,*result=NULL;
  char found=0;
  p=list;
  q=list->next;
  
  while (p!=q && found==0)
  {
    if (q->element==element)
    {
      found=1;
      result=q;
    }
    q=q->next;
  }
  return result;
}
// this function prints the list
void printlist(void)
{
  item *p;
  p=list->next;
  while (p!=list)
  {
    printf("%d ",p->element);
    p=p->next;
  }
}
void main(void)
{
  // creating first element of the list
  list=(item*)malloc(sizeof(item));
  list->element=0;
  list->next=list;
  list->prev=list;
  printf("Adding elements 1,2,4 & 5\n");
  Insert(list,5);
  Insert(list,4);
  Insert(list,2);
  Insert(list,1);
  printlist();
  printf("\nSearching element 2\n");
  item *tmp=Search(2);
  printf("Inserting element 3 after it\n");
  Insert(tmp,3);
  printlist();
  printf("\nDeleting element 1\n");
  Remove(list->next);
  printlist();
  // destroying list
  while (list->next!=list) Remove(list->next);
  free(list);
}

Binary tree
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
//////////////////////////////////////////////////////////////////////////////
//
//  Dynamic structures (binary tree)
//  (c) Johna Smith, 1996
//
//  Method description:
//               *
//            /     \
//           *       *
//         /   \   /   \
//        *     * *     *
//
//   From current element X left element is less than X and right element
// is greater than X. All elements in the must be different.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <alloc.h>
#include <conio.h>
#include <math.h>
struct item
{
  int element;
  item *left;
  item *right;
};
item *tree; // base element of the list
// this function searches element in the tree and returns 0 if element wasn't
// found and 1 - if element was found, result is address of element
char Search(int element, item** result)
{
  item *p,*q;
  char found=0;
  p=tree;
  if (tree!=NULL)
  do
  {
    q=p;
    if (p->element==element) found=1;
    else
    {
      q=p;
      if (element<p->element) p=p->left;
      else p=p->right;
    }
  }
  while (!found && p!=NULL);
  *result=q;
  return found;
}
// this function adds an element to the tree
void Add(int element)
{
  item *r,*s;
  if (Search(element,&r)==0)
  {
    s=(item*)malloc(sizeof(item));
    s->element=element;
    s->left=NULL;
    s->right=NULL;
    if (tree==NULL) tree=s; // if tree is empty make s=top of the tree
    else
    {
      if (element<r->element) r->left=s;
      else r->right=s;
    }
  }
}
// this is auxulary function for Remove procedure
void Del(item **r, item **q)
{
  item *tmp;
  if ((*r)->right==NULL)
  {
    (*q)->element=(*r)->element;
    *q=*r;
    *r=(*r)->left;
  } 
  else Del(&((*r)->right),q);
}
// this function removes element with value 'element' from the tree
void Remove(int element, item **d)
{
  item *q;
  if (*d==NULL)
  printf("There is not element %d in the tree.\n",element);
  else
  if (element<(*d)->element) Remove(element, &((*d)->left)); else
  if (element>(*d)->element) Remove(element, &((*d)->right)); else
  {
    // element found
    q=*d;
    if (q->right==NULL) *d=q->left; else
    if (q->left==NULL) *d=q->right; else
    Del(&(q->left),&q);
    free(q);
  }
}
// this function prints the tree
void printtree(item *t, int offset=40, int depth=2)
{
  gotoxy(offset,depth);
  cprintf("%d",t->element);
  if (t->left!=NULL) printtree(t->left,offset-pow(2,6-depth),depth+1);
  if (t->right!=NULL) printtree(t->right,offset+pow(2,6-depth),depth+1);
}
void main(void)
{
  item *tmp;
  // creating tree
  Add(100);
  Add(20);
  Add(120);
  Add(15);
  Add(50);
  Add(130);
  Add(30);
  Add(55);
  Add(28);
  Add(35);
  Add(60);
  Add(33);
  // printing tree
  clrscr();
  printf("Press a key to delete element 50...\n");
  printtree(tree);
  getch();
  clrscr();
  Remove(50,&tree);
  printtree(tree);
  gotoxy(1,20);
  // searching
  cprintf("Element 20 is%s found",(Search(20,&tmp)?"":"n't"));
  printf("\nElement 25 is%s found\n",(Search(25,&tmp)?"":"n't"));
  // removing all elements
  Remove(100,&tree);
  Remove(20,&tree);
  Remove(120,&tree);
  Remove(15,&tree);
  Remove(35,&tree);
  Remove(130,&tree);
  Remove(30,&tree);
  Remove(55,&tree);
  Remove(28,&tree);
  Remove(33,&tree);
  Remove(60,&tree);
}

Stack
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
//////////////////////////////////////////////////////////////////////////////
//
//  Dynamic structures (stack)
//  (c) Johna Smith, 1996
//
//  Method description:
//    Stack is the queue with LIFO structure (Last In First Out)
//    There is only one accesible element in the stack - its top
//    Two operations defined for stack - pushing element to stack
//    and popping element feom stack
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <alloc.h>
struct item
{
  int element;
  item *next;
};
item *stack=NULL; // base element of the stack
// this function pushes an to stack
void Push(int element)
{
  item *q;
  
  q=(item*)malloc(sizeof(item));
  q->element=element;
  q->next=stack;
  stack=q;
}
// this function pops an element from stack
int Pop(void)
{
  int result;
  item *q;
  
  if (stack==NULL)
    printf("!POP error: stack is EMPTY");
  else
  {
    result=stack->element;
    q=stack;
    stack=stack->next;
    free(q);
  }
  return result;
}
// this function prints the stack
void printstack(void)
{
  item *p;
  p=stack;
  while (p!=NULL)
  {
    printf("%d ",p->element);
    p=p->next;
  }
}
void main(void)
{
  printstack();
  printf("Pushing elements 1,2,3\n");
  Push(1);
  Push(2);
  Push(3);
  printstack();
  printf("\nPopping 2 elements\n");
  printf("%d ",Pop());
  printf("%d\n",Pop());
  printstack();
  // destroying stack
  while (stack!=NULL) Pop();
}

List
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
//////////////////////////////////////////////////////////////////////////////
//
//  Dynamic structures (list)
//  (c) Johna Smith, 1996
//
//  Method description:
//    * -> * -> * -> nil
//    A    B    C
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <alloc.h>
struct item
{
  int element;
  item *next;
};
item *list; // base element of the list
item *p; // current pointer
// this function adds an element to the end of list
void Add(int element)
{
  p->next=(item*)malloc(sizeof(item));
  p=p->next;
  p->element=element;
  p->next=NULL;
}
// this function removes element from list
// after element pointed by i
void Remove(item* i)
{
  item *tmp;
  if (i->next!=NULL)
  {
    tmp=i->next;
    i->next=(i->next)->next;
    free(tmp);
  }
}
// this function inserts element after
// element pointed by i
void Insert(item *i, int element)
{
  item *tmp;
  tmp=(item*)malloc(sizeof(item));
  tmp->element=element;
  tmp->next=i->next;
  i->next=tmp;
}
// this function searches element in the list
item* Search(int element)
{
  item *i,*result=NULL;
  char found=0;
  i=list->next;
  while (i!=NULL && found==0)
  {
    if (i->element==element)
    {
      found=1;
      result=i;
    }
    i=i->next;
  }
  return result;
}
// this function prints the list
void printlist(void)
{
  item *p;
  p=list->next;
  while (p!=NULL)
  {
    printf("%d ",p->element);
    p=p->next;
  }
}
void main(void)
{
  // creating first element of the list
  list=(item*)malloc(sizeof(item));
  list->element=0;
  list->next=NULL;
  p=list;
  printf("Adding elements 1,2 & 4\n");
  Add(1);
  Add(2);
  Add(4);
  printlist();
  printf("\nSearching element 2\n");
  item *tmp=Search(2);
  printf("Inserting element 3 after it\n");
  Insert(tmp,3);
  printlist();
  printf("\nDeleting element 1\n");
  Remove(list);
  printlist();
  // destroying list
  while (list->next!=NULL) Remove(list);
  free(list);
}


Сортировка массивов

Binary insertions
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
//////////////////////////////////////////////////////////////////////////////
//
//  Binary insertion (array sorting)
//  (c) Johna Smith, 1996
//
//  Method description:
//    Split our array into two parts - sorted (left) and unsorted (right)
//    Initially sorted part contains only one element - first array element
//    Take an element from the unsorted part and put it in the right place in
//    the sorted part (using binary search to find this right place). So
//    sorted part contains now two elements. Repeat this process until
//    no elements left in the unsorted part.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N=     10;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
int element; // value of the element we want to insert
unsigned int right,left,middle; // auxulary variables for binary search
void show_array(void) // displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void main(void)
{
  // Display unsorted array
  printf("Unsorted array: ");
  show_array();
  // Sorting
  for (int i=1;i<N;i++) // main loop
  {
    // Searching place for element i (binary search)
    left=0;
    right=i;
    element=array[i];
    while (left<right)
    {
      middle=(left+right)/2;
      if (array[middle]<=element) left=middle+1; else right=middle;
    }
    // Inserting element i into the right place
    for (int j=i;j>right;j--) array[j]=array[j-1];
    array[right]=element;
  }
  // Display sorted array
  printf("\nSorted array: ");
  show_array();
}
Straight insertion
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
//////////////////////////////////////////////////////////////////////////////
//
//  Straight insertion (array sorting)
//  (c) Johna Smith, 1996
//
//  Method description:
//    Split our array into two parts - sorted (left) and unsorted (right)
//    Initially sorted part contains only one element - first array element
//    Take an element from the unsorted part and put it in the right place in
//    the sorted part. So sorted part contains now two elements. Repeat
//    this process until no elements left in the unsorted part.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N = 10;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
int element; // value of the element we want to insert
unsigned int index; // index of the right position of the element
void show_array(void) // displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void main(void)
{
  // Display unsorted array
  printf("Unsorted array: ");
  show_array();
  // Sorting
  for (int i=1;i<N;i++) // main loop
  {
    index=0;
    // Searching place for element i
    while (array[index]<array[i]) index++;
    // Inserting element i into the right place
    element=array[i];
    for (int j=i;j>index;j--) array[j]=array[j-1];
    array[index]=element;
  }
  // Display sorted array
  printf("\nSorted array: ");
  show_array();
}

Heap sort
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
//////////////////////////////////////////////////////////////////////////////
//
//  Heap sort
//  (c) Johna Smith, 1996
//
//  Method description:
//    This algorithm sorts array using pyramids. The pyramid is the
//    sequence of characters h(L), h(L+1),...,h(R) with the following
//    properties:
//          h(i)<=h(2i);  h(i)<=h(2i+1), i=L..R/2
//                        h1
//                     /     \               This is an example of
//                   h2       h3                   pyramid
//                 /    \   /    \
//               h4     h5 h6     h7
//              .....................
//    There are two steps in this algorythm:
//    1) Building pyramid from the given array using Floyd method of
//    building pyramid "on the same place".
//    2) Sorting built pyramid:
//       Swap last pyramid element (x) with the top one and shift x
//       to the right place using Floyd method
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N = 10;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
int swp; // auxulary variable for swapping
int left,right; // indexes for left and right bounds of the sorted part of the array
void show_array(void) // this function displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
// this function shifts element in the pyramid to the right place
// it helps to build pyramid "on the same place" (using Floyd method)
void shift(int left, int right)
{
  int i,j;
  int swp;
  int element;
  i=left;
  j=2*left+1;
  element=array[left];
  if (j<right && array[j]<array[j+1]) j++;
  while (j<=right && element<array[j])
  {
    // swapping elements
    swp=array[i];
    array[i]=array[j];
    array[j]=swp;
    // now i-th element is on j-th place:
    i=j;
    element=array[i];
    // and j-th element is lower (see pyramid picture)
    j=2*j+1;
    if (j<right && array[j]<array[j+1]) j++;
  }
}
void main(void)
{
   // Displaying unsorted array
   printf("Unsorted array: ");
   show_array();
   // Sorting
   left=N/2;
   right=N-1;
   // building pyramid from array
   while (left>0)
   {
     left--;
     shift(left,N-1);
   }
   // sorting built pyramid
   while (right>0)
   {
     // swapping elements
     swp=array[0];
     array[0]=array[right];
     array[right]=swp;
     // shifting element to the right place
     right--;
     shift(0,right);
   }
   // Displaying sorted array
   printf("\nSorted array: ");
   show_array();
}

Straight selections
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Straight selections (array sorting)
//  (б) Johna Smith, 1996
//
//  Method description:
//    searching for smallest element in array, swapping with first element,
//    repeat this operation starting search from second element and so on
//    array become sorted when we'll start search from (n-1)-th element
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N = 10;  // number of elements in array
int array[N]={100,15,234,11,63,78,4,200,0,4};
int min; // smallest element in array
int index; // index of the smallest element in array
int swp; // auxulary variable for swapping
void show_array(void) // this function prints array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void main(void)
{
   // Printing unsorted array
   printf("Unsorted array: ");
   show_array();
   // Sorting
   for (int i=0;i<N-1;i++) // main loop
   {
     // searching for smallest element from i-th
     min=array[i]; // setting i-th elementas smallest
     index=i;
     for (int j=i+1;j<N;j++)
       if (array[j]<min)  // if there is element less than smallest
       {
          min=array[j];  // then remember its value
          index=j;       // and index
       }
     // swapping i-th and smallest element
     swp=array[index];
     array[index]=array[i];
     array[i]=swp;
   }
   // Printing sorted array
   printf("\nSorted array: ");
   show_array();
}

Shaker sort
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
//////////////////////////////////////////////////////////////////////////////
//
//  Shaker sort
//  (c) Johna Smith, 1996
//
//  Method description:
//    Advanced bubble sort.
//    There are two main steps:
//     1) Move greatest element to the end, moving forward and
//        remembering place of the last elements exchange
//        this place will be right bound of the unsorted part
//     2) Move smallest element to the beginning, moving backward and
//        remembering place of the last elements exchange
//        this place will be left bound of the unsorted part
//    Repeat this two steps while right bound is greater than left
//    This method is useful when almost all elements in the array
//    are sorted.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N=     10;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
int index; // auxulary variable (the place of last swap)
int swp; // auxulary variable for swapping
int left,right; // indexes for left and right bounds of the sorted part of the array
void show_array(void) // this function displays array
{
  for (int i=0;i<N;i++)
     printf("%d ",array[i]);
}
void main(void)
{
   // Displaying unsorted array
   printf("Unsorted array: ");
   show_array();
   // Sorting
   left=0;
   right=N-1;
   do
   {
     // moving smallest element to the beginning (moving backward)
     for (int i=right;i>left;i--)
     {
        if (array[i-1]>array[i])
        {
          // swapping a(i) and a(i+1)
          swp=array[i];
          array[i]=array[i-1];
          array[i-1]=swp;
          index=i; // save index
        }
     }
     left=index;
     // moving greatest element to the end (moving forward)
     for (i=left;i<right;i++)
     {
       if (array[i]>array[i+1])
       {
         // swapping a(i) and a(i+1)
         swp=array[i];
         array[i]=array[i+1];
         array[i+1]=swp;
         index=i; // save index
       }
     }
     right=index;
} while (left<right);
   // Displaying sorted array
   printf("\nSorted array: ");
   show_array();
}

Shell sort
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
//////////////////////////////////////////////////////////////////////////////
//
//  Shell sort (array sorting)
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Sort elements (0,4,8,12,...) (step=4) using straight insertions method
//    2) Sort elements (0,2,4,6,8,...) (step=2) using straight insertions method
//    3) Sort all elements (step=1) using straight insertions method
//    Shell sort is faster than simple straight insertions because at each
//    step more and more elements are already sorted (productivity increased
//    by ~300% for random array of 256 elements & about 700% (!) for random array
//    of 2048 elements (N. Wirth, Algoritms and data structures)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const T = 4;   // number of distances
const N = 20;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4,14,32,44,58,1,2,3,9,8,7}; // array of N integers
int interval[T]={9,5,3,1}; // last interval must be 1
int step; // current step
int element; // value of the element we want to insert
unsigned int index; // index of the right position of the element
void show_array(void) // displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void main(void)
{
  // Display unsorted array
  printf("Unsorted array: ");
  show_array();
  // Sorting
  for (int m=0;m<T;m++)  // main loop (by all distances)
  {
    step=interval[m];
    // sorting (0,step,2*step,3*step,...) elements using straight insertions
    for (int i=0;i<N;i+=step)
    {
      index=0;
      // Searching place for element i
      while (array[index]<array[i]) index+=step;
      // Inserting element i into the right place
      element=array[i];
      for (int j=i;j>index;j-=step) array[j]=array[j-step];
      array[index]=element;
    }
  }
  // Display sorted array
  printf("\nSorted array: ");
  show_array();
}

Quick sort (recursive)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
//////////////////////////////////////////////////////////////////////////////
//
//  Quick sort (recursive)
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Split array into two parts and remember middle element
//    2) Scan left part for element greater than middle
//    3) Scan right part for element less than middle
//    4) Swap these elements
//    So we have an array where all left elements are less than right elements
//    Apply these four steps to each part (left and right) of the array
//    until we have parts that contain only one element.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N=     10;  // array size
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
void show_array(void) // this function displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void sort(int left,int right)
{
  int i,j;
  int element; // auxulary variable for middle element in the interval
  int swp; // auxulary variable for swapping
  i=left;  // index for left part
  j=right;  // index for right part
  element=array[(left+right)/2]; // middle element
  do
  {
    while (array[i]<element) i++; // scanning left part
    while (element<array[j]) j--; // scanning right part
    if (i<=j)
    {
      // swapping elements
      swp=array[i];
      array[i]=array[j];
      array[j]=swp;
      i++; j--;
    }
  } while (i<=j);
  if (left<j) sort(left,j); // applying the same procedure to the left part
  if (i<right) sort(i,right); // applying the same procedure to the right part
}
void main(void)
{
  // Displaying unsorted array
  printf("Unsorted array: ");
  show_array();
  // Sorting
  sort(0,N-1);
  // Displaying sorted array
  printf("\nSorted array: ");
  show_array();
}

Quick sort (non-recursive)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
//////////////////////////////////////////////////////////////////////////////
//
//  Quick sort (non-recursive)
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Split array into two parts and remember middle element
//    2) Scan left part for element greater than middle
//    3) Scan right part for element less than middle
//    4) Swap these elements
//    So we have an array where all left elements are less than right elements
//    Apply these four steps to each part (left and right) of the array
//    until we have parts that contain only one element.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N = 10;         // array size
const STACK_SIZE = 4; // stack size must be >=log N
int array[N]={100,15,234,11,63,78,4,200,0,4}; // array of N integers
void show_array(void) // this function displays array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
struct {int left; int right;} stack[STACK_SIZE]; // stack emulation
unsigned int ss; // current stack size
int i,j;
int element; // auxulary variable for middle element in the interval
int swp; // auxulary variable for swapping
int left,right; // interval bounds
void main(void)
{
  // Displaying unsorted array
  printf("Unsorted array: ");
  show_array();
  // Sorting
  ss=1;
  stack[0].left=0;
  stack[0].right=N-1;
  do
  {
    // pushing from stack last request
    left=stack[ss-1].left;
    right=stack[ss-1].right;
    ss--;
    do
    {
      // splitting array[left]..array[right] interval
      i=left;  // index for left part
      j=right;  // index for right part
      element=array[(left+right)/2]; // middle element
      do
      {
        while (array[i]<element) i++; // scanning left part
        while (element<array[j]) j--; // scanning right part
        if (i<=j)
        {
          // swapping elements
          swp=array[i];
          array[i]=array[j];
          array[j]=swp;
          i++; j--;
        }
      } while (i<=j);
      // pushing request into stack
      // (we select longest part and push request for it to reduce stack size)
      if (j-left<right-i)
      {
        if (i<right) // push request for the right part (because it's longer than left)
        {
          ss++;
          stack[ss-1].left=i;
          stack[ss-1].right=right;
        }
        right=j; // continue sorting left part
      } else
      {
        if (left<j) // push request for the left part (because it's longer than right)
        {
          ss++;
          stack[ss-1].left=left;
          stack[ss-1].right=j;
        }
        left=i; // continue sorting right part
      }
    } while(left<right);
  } while(ss!=0);
  // Displaying sorted array
  printf("\nSorted array: ");
  show_array();
}

Exchanges
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
//////////////////////////////////////////////////////////////////////////////
//
//  Exchanges (array sorting)
//  (б) Johna Smith, 1996
//
//  Method description:
//    Using linear search find smallest i with the following property:
//    a(i)>a(i+1), swap a(i) and a(i+1) and repeat this process searching
//    from a(i+1). After one pass greatest number will be placed at right
//    place. We'll decrease amount of acting elements on each pass.
//    When 2 elements will left we'll say that array is sorted
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
const N=     10;  // number of array elements
int array[N]={100,15,234,11,63,78,4,200,0,4};
int swp; // auxulary variable for swapping
int m; // number of active elements on current pass
void show_array(void) // this function prints array
{
  for (int i=0;i<N;i++)
    printf("%d ",array[i]);
}
void main(void)
{
   // Printing unsorted array
   printf("Unsorted array: ");
   show_array();
   // Sorting
   m=N; // All elements acting at first pass
   while (m>1) // repeat passes while there is more than one element
   {
     for (int i=0;i<m-1;i++) // main loop
     {
       // searching for smallest element
       if (array[i]>array[i+1]) // if a(i)>a(i+1)
       {
         // swapping a(i) and a(i+1)
         swp=array[i];
         array[i]=array[i+1];
         array[i+1]=swp;
       }
     }
     m--; // decreasing amount of acting elements
   }
   // Printing sorted array
   printf("\nSorted array: ");
   show_array();
}


Сортировка файлов

Straight merge
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
//////////////////////////////////////////////////////////////////////////////
//
//  File sort (straight merge)
//  (c) Johna Smith, 1996
//
//  Method description:
//   This program reads data from file 1.dat and writes sorted data to 1.srt
//   Format of 1.dat: i1 i2 i3 i4 i5 i6 ..., where i(k) is integer
//   It also creates four temporary files: 1.tmp, 2.tmp, 3.tmp, 4.tmp
//   We use four files to sort given sequence of data.
//     1) Split file into 2 parts (A and B)
//     2) Combine A and B to C, sorting 2-element groups
//     3) Repeat these steps for file C, sorting 4-element groups,
//        8-elements,... while size of group less than file size.
//     We can optimize this method by uniting splitting and combining:
//     we just need to change destination on each step: combine one
//     group to C another to D, next group again to C and so on
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
int copy(FILE *A, FILE *B, FILE *C, int n)
{
  int a,b,q,r,m=0;
  q=r=0;
  // combining to C
  if (!feof(A)) fread(&a,sizeof(int),1,A);
  if (!feof(B)) fread(&b,sizeof(int),1,B);
  while (!feof(A) && !feof(B) && q<n && r<n)
  {
    if (a<b)
    {
      fwrite(&a,sizeof(int),1,C);m++;q++;
      if (q<n) fread(&a,sizeof(int),1,A);
    } 
    else
    {
      fwrite(&b,sizeof(int),1,C);m++;r++;
      if (r<n) fread(&b,sizeof(int),1,B);
    }
  }
  // writing remainders
  while (!feof(A) && q<n)
  {
    fwrite(&a,sizeof(int),1,C);m++;q++;
    if (q<n) fread(&a,sizeof(int),1,A);
  }
  while (!feof(B) && r<n)
  {
    fwrite(&b,sizeof(int),1,C);m++;r++;
    if (r<n) fread(&b,sizeof(int),1,B);
  }
  return m;
}
void main(void)
{
  FILE *A,*B,*C,*D,*input,*output;
  int N=0,c,m,p,l;
  // opening needed files (file1.dat must exist)
  input=fopen("1.dat","rb");
  output=fopen("1.srt","wb");
  l=0;
  A=fopen("1.tmp","w+b");
  B=fopen("2.tmp","w+b");
  C=fopen("3.tmp","w+b");
  D=fopen("4.tmp","w+b");
  // reading from source and splitting data into A & B
  while (!feof(input))
  {
    fscanf(input,"%d ",&c);
    fwrite(&c,sizeof(int),1,(N%2==0?A:B));
    N++;
  }
  rewind(A);
  rewind(B);
  p=1;
  while (p<N)
  {
    m=N;
    while (m!=0)
    {
      // uniting splitting and combining
      // here we combine from A and B to C and D
      m-=copy(A,B,C,p);
      m-=copy(A,B,D,p);
    }
    // swap A,B and C,D
    // now we'll combine from C and D and split to A and B
    fclose(A);
    fclose(B);
    fclose(C);
    fclose(D);
    A=fopen((l?"1.tmp":"3.tmp"),"rb");
    B=fopen((l?"2.tmp":"4.tmp"),"rb");
    C=fopen((!l?"1.tmp":"3.tmp"),"wb");
    D=fopen((!l?"2.tmp":"4.tmp"),"wb");
    l=!l;
    // groups will be larger twice
    p*=2;
  }
  fread(&c,sizeof(int),1,A);
  while (!feof(A))
  {
    fprintf(output,"%d ",c);
    fread(&c,sizeof(int),1,A);
  }
  fclose(A);
  fclose(B);
  fclose(C);
  fclose(D);
  fclose(input);
  fclose(output);
}

Natural merge
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
//////////////////////////////////////////////////////////////////////////////
//
//  File sort (natural merge)
//  (c) Johna Smith, 1996
//
//  Method description:
//   This program reads data from file 1.dat and writes sorted data to 1.srt
//   Format of 1.dat: i1 i2 i3 i4 i5 i6 ..., where i(k) is integer
//   It also creates four temporary files: 1.tmp, 2.tmp, 3.tmp
//   We use three files to sort given sequence of data.
//
//   Series - is a sequence of numbers with the following property: a(i)<=a(i+1)
//   1) Divide copy of the given file C into parts A and B, writing first
//      series from C to A, second one to B, third - to C, and so on
//   2) Combine files A and B into C, sorting distributed series
//   3) Repeat first two steps until there is more than one series
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
void main(void)
{
  FILE *A,*B,*C,*input,*output;
  int sorted,a,b,c,old;
  // opening needed files (file1.dat must exist)
  input=fopen("1.dat","rb");
  output=fopen("1.srt","wb");
  A=fopen("1.tmp","w+b");
  B=fopen("2.tmp","w+b");
  C=fopen("3.tmp","w+b");
  // reading from source and copying data into C
  while (!feof(input))
  {
    fscanf(input,"%d ",&c);
    fwrite(&c,sizeof(int),1,C);
  }
  rewind(C);
  do
  {
    // reopeneing files A and B
    fclose(A);
    fclose(B);
    A=fopen("1.tmp","w+b");
    B=fopen("2.tmp","w+b");
    rewind(C);
    // distributing series from C to A and B
    fread(&c,sizeof(int),1,C);
    sorted=1;
    while (!feof(C))
    {
      // series to A
      do
      {
         old=c;
         fwrite(&c,sizeof(int),1,A);
         fread(&c,sizeof(int),1,C);
      } while (!feof(C) && old<=c);
      old=c;
      // series to B
      while (!feof(C) && old<=c)
      {
        old=c;
        fwrite(&c,sizeof(int),1,B);
        fread(&c,sizeof(int),1,C);
        sorted=0; // if there are more than 1 series then sequence isn't sorted
      }
    }
    // reopening file C
    fclose(C);
    C=fopen("3.tmp","w+b");
    rewind(A);
    rewind(B);
    // combining A and B back to C
    fread(&a,sizeof(int),1,A);
    fread(&b,sizeof(int),1,B);
    while (!feof(A) && !feof(B))
    {
      if (a<b)
      {
        old=a;
        fwrite(&a,sizeof(int),1,C);
        fread(&a,sizeof(int),1,A);
        if (feof(A) || a<old)
        do
        {
          old=b;
          fwrite(&b,sizeof(int),1,C);
          fread(&b,sizeof(int),1,B);
        } while(!feof(B) && b>=old);
      } 
      else
      {
        old=b;
        fwrite(&b,sizeof(int),1,C);
        fread(&b,sizeof(int),1,B);
        if (feof(B) || b<old)
        do
        {
          old=a;
          fwrite(&a,sizeof(int),1,C);
          fread(&a,sizeof(int),1,A);
        } while(!feof(A) && a>=old);
      }
    }
    // copying remainder from A (or B) to C
    while (!feof(A))
    {
      fwrite(&a,sizeof(int),1,C);
      fread(&a,sizeof(int),1,A);
    }
    while (!feof(B))
    {
      fwrite(&b,sizeof(int),1,C);
      fread(&b,sizeof(int),1,B);
    }
  } while (!sorted);  // sort until there will be only one series
  // writing sorted sequence
  rewind(C);
  fread(&c,sizeof(int),1,C);
  while (!feof(C))
  {
    fprintf(output,"%d ",c);
    fread(&c,sizeof(int),1,C);
  }
  // closing all opened files
  fclose(A);
  fclose(B);
  fclose(C);
  fclose(input);
  fclose(output);
}

Balanced merge
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
//////////////////////////////////////////////////////////////////////////////
//
//  File sort (balanced merge)
//  (c) Johna Smith, 1996
//
//  Method description:
//   This program reads data from file 1.dat and writes sorted data to 1.srt
//   Format of 1.dat: i1 i2 i3 i4 i5 i6 ..., where i(k) is integer
//   It also creates four temporary files: 1.tmp, 2.tmp, 3.tmp
//   We use N files to sort given sequence of data.
//
//   Series - is a sequence of numbers with the following property: a(i)<=a(i+1)
//   1) Divide copy of the given file C into N files, writing first
//      series from C to F1, second one to F2, third - to F3, and so on
//   2) Combine files F1,F2,F3,...,Fn, sorting distributed series
//   3) Repeat first two steps until there is more than one series
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#define N   10  // number of temporary files (must be even and greater than 2)
#define Nh  N/2 // must be an integer
struct ExtFile
{
  FILE *F;
  int first; // current element of this file
  char eor; // indicates end of series
};
// this function copies data element from x to y, updating (.first) members
void copy(ExtFile *x, ExtFile *y)
{
  y->first=x->first;
  fwrite(&y->first,sizeof(int),1,y->F);
  fread(&x->first,sizeof(int),1,x->F);
  x->eor=x->first<y->first;
}
void main(void)
{
  FILE *input,*output;
  ExtFile F[N];
  int l,old,j,mx,tx,min,x,t[N],ta[N],k1,k2;
  char name[13];
  // opening needed files (file1.dat must exist)
  input=fopen("1.dat","rb");
  output=fopen("1.srt","wb");
  for (int i=0;i<N;i++)
  {
    gcvt(i,8,name);
    strcat(name,".tmp");
    F[i].F=fopen(name,"w+b");
    F[i].eor=0;
  }
  // distributing given data
  j=Nh;
  l=0;
  fscanf(input,"%d ",&x);
  while (!feof(input))
  {
    if (j<Nh-1) j++; else j=0;
    do
    {
      old=x;
      fwrite(&x,sizeof(int),1,F[j].F);
      fscanf(input,"%d ",&x);
    } while (!feof(input) && x>=old);
    l++;
  }
  if (x<old)
  {
    if (j<Nh-1) j++; else j=0;
    l++;
  }
  fwrite(&x,sizeof(int),1,F[j].F);
  for (i=0;i<N;i++) t[i]=i;
  do
  {
    // merging t[0]..t[Nh-1] to t[Nh]..t[N-1]
    if(l<Nh) k1=l-1; else k1=Nh-1;
    for(i=0;i<=k1;i++)
    {
      rewind(F[t[i]].F);
      ta[i]=t[i];
      fread(&F[t[i]].first,sizeof(int),1,F[t[i]].F);
    }
    l=0; // number of input series
    j=Nh; // index of input sequence
    // merging input series to t[j]
    do
    {
      l++;
      k2=k1;
      // selecting minimal element from F1..Fn
      do
      {
        i=1;
        mx=0;
        min=F[ta[0]].first;
        while (i<=k2)
        {
          x=F[ta[i]].first;
          if(x<min)
          {
            min=x;
            mx=i;
          }
          i++;
        }
        copy(&F[ta[mx]],&F[t[j]]);
        if (feof(F[ta[mx]].F))
        {
          // excluding file
          fclose(F[ta[mx]].F);
          gcvt(ta[mx],8,name);
          strcat(name,".tmp");
          F[ta[mx]].F=fopen(name,"w+b");
          F[ta[mx]].eor=0;
          ta[mx]=ta[k2];
          ta[k2]=ta[k1];
          k1--;
          k2--;
        } else
        if (F[ta[mx]].eor)
        {
          // closing series
          tx=ta[mx];
          ta[mx]=ta[k2];
          ta[k2]=tx;
          k2--;
        }
      } while(k2>=0);
      if(j<N-1) j++; else j=Nh;
    } while(k1>=0);
    for(i=0;i<Nh;i++)
    {
      tx=t[i];
      t[i]=t[i+Nh];
      t[i+Nh]=tx;
    }
  } while (l!=1);
  // writing sorted sequence
  rewind(F[0].F);
  fread(&x,sizeof(int),1,F[0].F);
  while (!feof(F[0].F))
  {
    fprintf(output,"%d ",x);
    fread(&x,sizeof(int),1,F[0].F);
  }
  // closing all opened files
  for (i=0;i<N;i++) fclose(F[i].F);
  fclose(input);
  fclose(output);
}

Polyphase merge
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
//////////////////////////////////////////////////////////////////////////////
//
//  File sort (polyphase merge)
//  (c) Johna Smith, 1996
//
//  Method description:
//   This program reads data from file 1.dat and writes sorted data to 1.srt
//   Format of 1.dat: i1 i2 i3 i4 i5 i6 ..., where i(k) is integer
//   It also creates four temporary files: 1.tmp, 2.tmp, 3.tmp
//   We use N files to sort given sequence of data.
//
//   Series - is a sequence of numbers with the following property: a(i)<=a(i+1)
//   1) Divide copy of the given file C into N-1 files, writing first
//      series from C to F1, second one to F2, third - to F3, and so on
//   2) Combine series from files F1,F2,F3,...,F(n-1) to Fn,
//      sorting distributed series
//   3) When we reach the end of one of F1..F(n-1) files we turn sequences
//      and will combine series to this new empty file. For example if we reach
//      end of F3 we will gather data from F1,F2,F4,..Fn to F3
//   4) Repeat first three steps until there is more than one series
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#define N   6  // number of temporary files (must be greater than 2)
struct ExtFile
{
  FILE *F;
  int first; // current element of this file
  char eor; // indicates end of series
};
  FILE *output;
  ExtFile input;
  ExtFile F[N];  // array of files
  int level,i,j,k,mx,tx,xmin,x,t[N],ta[N],a[N],d[N],z,tn,dn;
  char name[13];
// this functions sets new value of variable j
// on the next step of data distribution we'll write to sequence #j
void select(void)
{
  if (d[j]<d[j+1]) j++;  // d[j] shows number of series in sequence #j
  else
  {
    if (d[j]==0)
    {
      level++;
      z=a[0];
      for (int i=0;i<N-1;i++)
      {
        d[i]=z+a[i+1]-a[i];
        a[i]=z+a[i+1];
      }
    }
    j=0;
  }
  d[j]--;
}
// this function copies data element from x to y, updating (.first) members
void copy(ExtFile *x, ExtFile *y)
{
  y->first=x->first;
  fwrite(&y->first,sizeof(int),1,y->F);
  fread(&x->first,sizeof(int),1,x->F);
  x->eor=(x->first<y->first || feof(x->F));
}
// this function copies series from input to f[j]
void copyrun(void)
{
  do
  {
    F[j].first=input.first;
    fwrite(&F[j].first,sizeof(int),1,F[j].F);
    fscanf(input.F,"%d",&input.first);
    input.eor=input.first<F[j].first;
  } while (!input.eor && !feof(input.F));
  if (feof(input.F))
  {
    if (input.first<F[j].first) select();
    F[j].first=input.first;
    fwrite(&F[j].first,sizeof(int),1,F[j].F);
    d[j]++;
  }
}
void main(void)
{
  // opening required files (file 1.dat must exist)
  input.F=fopen("1.dat","rb");
  fscanf(input.F,"%d ",&input.first);
  F[t[N-1]].eor=0;
  output=fopen("1.srt","wb");
  for (int i=0;i<N-1;i++)
  {
    gcvt(i,8,name);
    strcat(name,".tmp");
    F[i].F=fopen(name,"w+b");
    fread(&F[i].first,sizeof(int),1,F[i].F);
    F[i].eor=0;
    a[i]=1;
    d[i]=1;
  }
  // initializing variables
  level=1;
  j=0;
  a[N-1]=0;
  d[N-1]=0;
  // distributing data from input to first N-1 sequences
  do
  {
    select();
    copyrun();
  } while (!feof(input.F) && j!=N-2);
  while (!feof(input.F))
  {
    select();
    if (F[j].first<=input.first)
    {
      copyrun();
      if (feof(input.F)) d[j]++; else copyrun();
    } else copyrun();
  }
  // preparing first N-1 series for reading
  for (i=0;i<N-1;i++)
  {
    t[i]=i;
    rewind(F[i].F);
    fread(&F[i].first,sizeof(int),1,F[i].F);
    F[i].eor=0;
  }
  t[N-1]=N-1;
  // main loop
  do
  {
    // gathering data from t[0]..t[N-2] into t[N-1]
    z=a[N-2];
    d[N-1]=0;
    // preparing file to write
    fclose(F[t[N-1]].F);
    gcvt(t[N-1],8,name);
    strcat(name,".tmp");
    F[t[N-1]].F=fopen(name,"w+b");
    F[t[N-1]].eor=0;
    do
    {
      // gathering one series
      k=-1;
      for (i=0;i<N-1;i++)
      {
        if (d[i]>0) d[i]--;
        else
        {
          k++;
          ta[k]=t[i];
        }
      }
      if (k==-1) d[N-1]++;
      else
      {
        // gathering one real series from t[0]..t[k] into t[N-1]
        do
        {
          i=0;
          mx=0;
          xmin=F[ta[i]].first;
          // selecting series with minimal element to place it first
          while (i<k)
          {
            i++;
            x=F[ta[i]].first;
            if (x<xmin)
            {
              xmin=x;
              mx=i;
            }
          }
          // copying smallest element to output series
          copy(&F[ta[mx]],&F[t[N-1]]);
          // if we reach the end of the current input sequence
          // or end of series in it then close this source
          if (F[ta[mx]].eor)
          {
            // closing this source
            ta[mx]=ta[k];
            k--;
          }
        } while (k!=-1);
      }  
      z--;
    } while (z!=0);
    // swapping sequences
    rewind(F[t[N-1]].F);
    fread(&F[t[N-1]].first,sizeof(int),1,F[t[N-1]].F);
    F[t[N-1]].eor=0;
    tn=t[N-1];
    dn=d[N-1];
    z=a[N-2];
    for(i=N-1;i>0;i--)
    {
      t[i]=t[i-1];
      d[i]=d[i-1];
      a[i]=a[i-1]-z;
    }
    t[0]=tn;
    d[0]=dn;
    a[0]=z;
    // preparing file F[t[N-1]] to write
    fclose(F[t[N-1]].F);
    gcvt(t[N-1],8,name);
    strcat(name,".tmp");
    F[t[N-1]].F=fopen(name,"w+b");
    F[t[N-1]].eor=0;
    level--;
  } while (level!=0); // while there are more than one series
  // writing sorted sequence from F[t[0]]
  rewind(F[t[0]].F);
  fread(&x,sizeof(int),1,F[t[0]].F);
  while (!feof(F[t[0]].F))
  {
    fprintf(output,"%d ",x);
    fread(&x,sizeof(int),1,F[t[0]].F);
  }
  // closing all opened files
  for (i=0;i<N;i++) fclose(F[i].F);
  fclose(input.F);
  fclose(output);
}


Методы поиска

Binary search in sorted string
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
//////////////////////////////////////////////////////////////////////////////
//
//  Binary search in sorted string
//  (c) Johna Smith, 1996
//
//  Method description:
//     1) Split string into two parts
//     2) If the middle symbol of the string is greater than searched one -
//        select left part else select right part
//     3) Split selected part into two parts (see 1)
//   When in the selected part will be only one symbol compare it with
//   searched one.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
char string[]="0245789AAABDFHJKLXYZ"; // sorted ASCIIZ string
char symbol='J'; // symbol to find
unsigned int left,right; // bounds of the area where we search
unsigned int middle; // center of the area where we search
unsigned int length=sizeof(string)/sizeof(char); // string length
void main(void)
{
  printf("String: %s\n",string);
  left=0;
  right=length-1;
  // Here we avoid checking condition 'string[middle]==symbol'
  // because maximum number of steps is log N (where N is length
  // of the string) and some auxulary steps will be more effective
  // than comparing elements each step.
  while (left<right)
  {
    middle=(left+right)/2;
    if (string[middle]<symbol) left=middle+1; else right=middle;
  }
  // We'll stop when left==right => if there is searched symbol in
  // the string both bounds will point at it
  if (string[right]==symbol) printf("Symbol '%c' found at position %d.",symbol,right+1);
  else printf("Symbol '%c' not found.",symbol);
}

Substring search. Boyer-Moore algorithm
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
//////////////////////////////////////////////////////////////////////////////
//
//  Substring search. Boyer-Moore algorithm
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Build auxulary array D. D contains 256 elements (it's
//    amount of symbols in ASCII table. D(i) equals to length
//    of the substring if it doesn't contain character with code i.
//    And if substring contains character with code i then D(i)
//    equals to the distance between end of string and position of
//    this character closest to the end of string.
//    2) Scan string from the start, comparing it with the substring
//    If diferrence found we can skip next d(string(i)) positions,
//    where string(i) is the code of i-th character of the string and
//    i is the current position in the string.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
char string[]="Unsorted string. ABCBBCBCA."; // ASCIIZ string
char substring[]="string"; // substring to find
const unsigned int string_length=sizeof(string)/sizeof(char)-1; // string length
const unsigned int substring_length=sizeof(substring)/sizeof(char)-1; // substring length
int i,j,k;
int d[256]; // auxulary array
void main(void)
{
  printf("String: %s\n",string);
  // Analysing substring, building array d
  for (int i=0;i<256;i++) d[i]=substring_length;
  for (i=0;i<substring_length-1;i++) d[substring[i]]=substring_length-i-1;
  // Searching
  i=substring_length;
  do
  {
    j=substring_length; k=i;
    do // search difference between string and substring
    {
      k--;
      j--;
    } while (j>=0 && string[k]==substring[j]);
    i+=d[string[i-1]]; // skip next d(string(i)) positions
  } while (j>=0 && i<string_length);
  // condition j<0 means that substring is found
  if (j<0) printf("Substring \"%s\" found.",substring);
  else printf("Substring \"%s\" not found.",substring);
}

Linear substring search
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
//////////////////////////////////////////////////////////////////////////////
//
//  Linear substring search
//  (c) Johna Smith, 1996
//
//  Method description:
//    Simple search method.
//    1) Compare string and substring from first symbol of the string
//    2) If difference found compare string and substring from
//       second symbol of the string
//    Repeat these steps until no differencies found or end of string reached
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
char string[]="Unsorted string"; // ASCIIZ string
char substring[]="ring"; // substring to find
unsigned int position; // comparing from this position of the string
unsigned int subposition; // compared position of the substring
unsigned int string_length=sizeof(string)/sizeof(char); // string length
unsigned int substring_length=sizeof(substring)/sizeof(char); // substring length
void main(void)
{
  printf("String: %s\n",string);
  position=0;
  // finish compare when all symbols of substring compared (substring found)
  // or there left less than substring_length uncompared symbols
  // in the string (substring not found)
  while (subposition!=substring_length && position!=string_length-substring_length+1)
  {
    subposition=0; // comparing from the beginning of the substring
    while (subposition<substring_length &&
              string[position+subposition]==substring[subposition])
      subposition++;
    // compare from next position of the string
    position++;
  }
  if (subposition==substring_length) printf("Substring \"%s\" found from position %d.",substring,position);
  else printf("Substring \"%s\" not found.",substring);
}

Substring search. Knuth-Morris-Pratt algorithm
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Substring search. Knuth-Morris-Pratt algorithm
//  (c) Johna Smith, 1996
//
//  Method description:
//    1) Build auxulary array D. D contains "substring_length" elements.
//    D(i) equals to the length of longest sequence of characters
//    of the substring with such properties:
//         ss(0)..ss(d(i)-1)==ss(i-d(i))..p(i-1);  (ss is a substring)
//         p(d(i))!=p(i);
//    2) Scan string from the start, comparing it with the substring
//    If diferrence found we can skip next d(string(i)) positions,
//    where string(i) is the code of i-th character of the string and
//    i is the current position in the string.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
char string[]="Unsorted string. ABCBBCBCA."; // ASCIIZ string
char substring[]="BCB"; // substring to find
const unsigned int string_length=sizeof(string)/sizeof(char)-1; // string length
const unsigned int substring_length=sizeof(substring)/sizeof(char)-1; // substring length
int i,j,k;
int d[substring_length]; // auxulary array
void main(void)
{
  printf("String: %s\n",string);
  // Analysing substring, building array d
  j=0; k=-1; d[0]=-1;
  while (j<substring_length-1)
  {
    while (k>=0 && substring[j]!=substring[k]) k=d[k];
    j++; k++;
    if (substring[j]==substring[k]) d[j]=d[k]; else d[j]=k;
  }
  // Searching
  i=j=k=0;
  while (j<substring_length && i<string_length)
  {
    // here j is the current position in the substring and i is
    // the current position in the string
    while (j>=0 && string[i]!=substring[j]) j=d[j]; //skip next d[j] positions
    i++; j++;
  }
  // Condition j==substring length means that substring is found
  if (j==substring_length) printf("Substring \"%s\" found.",substring);
  else printf("Substring \"%s\" not found.",substring);
}

Linear search in the string
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Linear search
//  (c) Johna Smith, 1996
//
//  Method description:
//    Simple search method. Take first element - if it is what we
//   need - stop searh else take next element.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
char string[]="Unsorted string"; // ASCIIZ string
char symbol='d'; // symbol to find
char *p;
unsigned int position; // position of the found symbol
unsigned int length=sizeof(string)/sizeof(char); // string length
void main(void)
{
  printf("String: %s\n",string);
  // changing last (terminating) symbol of the string
  // to the symbol that we search allows us to avoid
  // checking condition 'end of string' - we always
  // will find symbol in the modified string
  string[length-1]=symbol;
  for (p=string,position=1;*p!=symbol;p++,position++);
  string[length-1]='\0'; //restoring last symbol
  
  if (position!=length) printf("Symbol '%c' found at position %d.",symbol,position);
  else printf("Symbol '%c' not found.",symbol);
}


Матрицы

Определитель
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating det A, where A is matrix N x N
//  (c) Johna Smith, 1996
//
//  Method description:
//    This method based on two statements:
//    1) det A = a11, if A is matrix 1 x 1
//    2) |a11 a12 a13 ... a1n|                       |a11 a12 .. a1(i-1) a1(i+1) .. a1n |
//       |a21 a22 a23 ... a2n|    n     i+1          |a21 a22 .. a2(i-1) a2(i+1) .. a2n |
// det A=|a31 a32 a33 ... a3n| = SUM (-1)  a1i * det |a31 a32 .. a3(i-1) a3(i+1) .. a3n |
//       |...................|   i=1                 |................................  |
//       |...................|                       |an1 an2 .. an(i-1) an(i+1) .. ann |
//       |an1 an2 an3 ... ann|
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <alloc.h>
#define N       5  // matrix size
int matrix[N][N]={{1,2,3,4,5},{2,1,3,4,5},{3,2,1,4,5},{4,3,2,1,5},{5,4,3,2,1}};
void show_matrix(int *matrix, int size) // this functions displays matrix
{
 for(int i=0;i<size;i++)
 {
   for(int j=0;j<size;j++)
     printf("%d ",*(matrix+i*size+j));
   printf("\n");
 }
 printf("\n");
}
int Det(int *matrix, int const size)
{
  int *submatrix;
  int det=0;
  if (size>1)
  {
    submatrix=(int *)malloc(sizeof(int)*(size-1)*(size-1));
    for (int i=0;i<size;i++)
    {
      // creating new array (submatrix)
      for (int j=0;j<size-1;j++)
        for (int k=0;k<size-1;k++)
          *(submatrix+j*(size-1)+k)=*(matrix+(j+1)*size+(k<i?k:k+1));
      // calling recursively function Det using submatrix as a parameter
      // and adding the result to final value
      det+=*(matrix+i)*(i/2.0==i/2?1:-1)*Det(submatrix,size-1);
    }
    free(submatrix);
  } else det=*matrix;
  return det;
}
void main(void)
{
  // show given matrix
  show_matrix(matrix[0],N);
  // calculating determinante and displaying the result
  printf("det A = %d",Det(matrix[0],N));
}

Обращение
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
//////////////////////////////////////////////////////////////////////////////
//
//  Matrix operations (inversion)
//  (c) Johna Smith, 1996
//
//  Given: A (N x N), det A != 0                 -1
//  This algorithm inverts matrix A and returns A  .
//      -1
//   A*A  = E
//  We simply convert matrix A into matrix E, and do the same operations
//  over matrix B (initially B=E). Finally A=E, B=A^(-1).
//
//////////////////////////////////////////////////////////////////////////////
#define N               3
#include <stdio.h>
float a[N][N]={{4,8,0},{8,8,8},{2,0,1}};
void Invert(float *matrix)
{
  float e[N][N];
  // initializing matrix e
  for (int i=0;i<N;i++)
    for (int j=0;j<N;j++)
      e[i][j]=(i==j?1:0);
  // converting matrix to e
  for(i=0;i<N;i++)
  {
    // normalizing row (making first element =1)
    float tmp=*(matrix+i*N+i);
    for(int j=N-1;j>=0;j--)
    {
      e[i][j]/=tmp;
      *(matrix+i*N+j)/=tmp;
    }
    // excluding i-th element from each row except i-th one
    for(j=0;j<N;j++)
      if (j!=i)
      {
        tmp=*(matrix+j*N+i);
        for(int k=N-1;k>=0;k--)
        {
          e[j][k]-=e[i][k]*tmp;
          *(matrix+j*N+k)-=*(matrix+i*N+k)*tmp;
        }
      }
  }
  // now e contains inverted matrix so we need only to copy e to matrix
  for(i=0;i<N;i++)
    for(int j=0;j<N;j++)
      *(matrix+i*N+j)=e[i][j];
}
void show_matrix(float *matrix) // this functions displays matrix
{
 for(int i=0;i<N;i++)
 {
   for(int j=0;j<N;j++)
     printf("%f ",*(matrix+i*N+j));
   printf("\n");
 }
 printf("\n");
}
void main(void)
{
  // display matrix A
  show_matrix(a[0]);
  // Invert it
  Invert(a[0]);
  // display the inverted matrix
  show_matrix(a[0]);
  // Invert it again
  Invert(a[0]);
  // display the inversion of the inverted matrix
  show_matrix(a[0]);
}

Сложение
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
//////////////////////////////////////////////////////////////////////////////
//
//  Matrix operations (addition)
//  (c) Johna Smith, 1996
//
//  Given: A,B - matrixes M x N
//  This algorithm sum these matrixes and display the result.
//
//////////////////////////////////////////////////////////////////////////////
#define M    5
#define N    6
#include <stdio.h>
int a[N][M]={{1,2,3,4,5},{5,4,3,2,1},{1,2,3,4,5},{5,4,3,2,1},{1,2,3,4,5},{0,1,0,1,0}};
int b[N][M]={{5,4,3,2,1},{1,2,3,4,5},{5,4,3,2,1},{1,2,3,4,5},{5,4,3,2,1},{1,0,1,0,1}};
int c[N][M];
void show_matrix(int *matrix) // this functions displays matrix
{
 for(int i=0;i<N;i++)
 {
    for(int j=0;j<M;j++)
      printf("%d ",*(matrix+i*M+j));
    printf("\n");
 }
 printf("\n");
}
void main(void)
{
  // display matrixes A and B
  show_matrix(a[0]);
  show_matrix(b[0]);
  // sum these matrixes
  for(int i=0;i<N;i++)
    for(int j=0;j<M;j++)
      c[i][j]=a[i][j]+b[i][j];
  // display the sum
  show_matrix(c[0]);
}

Умножение
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
//////////////////////////////////////////////////////////////////////////////
//
//  Matrix operations (multiplication)
//  (c) Johna Smith, 1996
//
//  Given: A (M x N) and B (N x P)
//  This algorithm multiplies these matrixes and display the result.
//
//////////////////////////////////////////////////////////////////////////////
#define M     3
#define N     4
#define P     2
#include <stdio.h>
int a[M][N]={{2,0,3,1},{5,1,2,0},{0,0,4,1}};
int b[N][P]={{1,3},{2,1},{4,0},{3,5}};
int c[M][P];
void show_matrix(int *matrix, int n, int m) // this functions displays matrix
{
 for(int i=0;i<n;i++)
 {
   for(int j=0;j<m;j++)
     printf("%d ",*(matrix+i*m+j));
   printf("\n");
 }
 printf("\n");
}
void main(void)
{
  // display matrixes A and B
  show_matrix(a[0],M,N);
  show_matrix(b[0],N,P);
  // substract these matrixes
  for(int i=0;i<M;i++)
    for(int j=0;j<P;j++)
    {
      c[i][j]=0;
      for (int k=0;k<N;k++) c[i][j]+=a[i][k]*b[k][j];
    }
  // display the difference
  show_matrix(c[0],M,P);
}
LK
Заблокирован
19.05.2013, 15:58  [ТС]     Коллекция алгоритмов от Johna Smith #3
Вычитание
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
//////////////////////////////////////////////////////////////////////////////
//
//  Matrix operations (substraction)
//  (c) Johna Smith, 1996
//
//  Given: A,B - matrixes M x N
//  This algorithm substract these matrixes and display the result.
//
//////////////////////////////////////////////////////////////////////////////
#define M   5
#define N   6
#include <stdio.h>
int a[N][M]={{4,4,4,4,4},{5,4,3,2,1},{1,2,3,4,5},{5,4,3,2,1},{1,2,3,4,5},{0,1,0,1,0}};
int b[N][M]={{3,4,4,4,4},{5,3,3,2,1},{1,2,2,4,5},{5,4,3,1,1},{1,2,3,4,4},{0,1,0,1,0}};
int c[N][M];
void show_matrix(int *matrix) // this functions displays matrix
{
 for(int i=0;i<N;i++)
 {
    for(int j=0;j<M;j++)
      printf("%d ",*(matrix+i*M+j));
    printf("\n");
 }
 printf("\n");
}
void main(void)
{
  // display matrixes A and B
  show_matrix(a[0]);
  show_matrix(b[0]);
  // substract these matrixes
  for(int i=0;i<N;i++)
    for(int j=0;j<M;j++)
      c[i][j]=a[i][j]-b[i][j];
  // display the difference
  show_matrix(c[0]);
}

Транспонирование
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Matrix operations (transponating)
//  (c) Johna Smith, 1996
//
//  Given: A (N x N),                                 T
//  This algorithm transponates matrix A and returns A  .
//          T
//   A(ji)=A (ij)
//
//////////////////////////////////////////////////////////////////////////////
#define N  3
#include <stdio.h>
float a[N][N]={{4,8,0},{8,8,8},{2,0,1}};
void Transponate(float *matrix)
{
  float swp;
  for (int i=0;i<N;i++)
    for (int j=i+1;j<N;j++)
    {
      swp=*(matrix+i*N+j);
      *(matrix+i*N+j)=*(matrix+j*N+i);
      *(matrix+j*N+i)=swp;
    }
}
void show_matrix(float *matrix) // this functions displays matrix
{
 for(int i=0;i<N;i++)
 {
   for(int j=0;j<N;j++)
     printf("%f ",*(matrix+i*N+j));
   printf("\n");
 }
 printf("\n");
}
void main(void)
{
  // display matrix A
  show_matrix(a[0]);
  // Transponate it
  Transponate(a[0]);
  // display transponated matrix
  show_matrix(a[0]);
  // Transponate it again
  Transponate(a[0]);
  // display the transponation of the transponated matrix
  show_matrix(a[0]);
}


Операции с комплексными величинами

Сложение и вычитание в алгебраической форме
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (addition and substraction in a+bi form)
//  (c) Johna Smith, 1996
//
//  Given: z1,z2 - complex values
//  z1=a1+i*b1, z2=a2+i*b2
//  z1+z2=(a1+a2)+i*(b1+b2)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Add(complex a, complex b)
{
  complex c;
  c.re=a.re+b.re;
  c.im=a.im+b.im;
  return c;
}
complex Substract(complex a, complex b)
{
  complex c;
  c.re=a.re-b.re;
  c.im=a.im-b.im;
  return c;
}
complex a={2,3},b={4,6};
void main(void)
{
  // addition
  show_complex(a);
  printf(" + ");
  show_complex(b);
  printf(" = ");
  show_complex(Add(a,b));
  // substraction
  printf("\n");
  show_complex(b);
  printf(" - ");
  show_complex(a);
  printf(" = ");
  show_complex(Substract(b,a));
}

Деление в алгебраической форме
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (division in a+bi form)
//  (c) Johna Smith, 1996
//
//  Given: z1,z2 - complex values
//  z1=a1+i*b1, z2=a2+i*b2
//  z1/z2=(a1*a2-b1*b2)/(a2^2+b2^2)+i*(b1*a2-a1*b2)/(a2^2+b2^2)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Div(complex a, complex b)
{
  complex c;
  c.re=(a.re*b.re+a.im*b.im)/(b.re*b.re+b.im*b.im);
  c.im=(a.im*b.re-a.re*b.im)/(b.re*b.re+b.im*b.im);
  return c;
}
complex a={2,3},b={-1,2};
void main(void)
{
  show_complex(a);
  printf(" / ");
  show_complex(b);
  printf(" = ");
  show_complex(Div(a,b));
}

Умножение в алгебраической форме
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (multiplication in a+bi form)
//  (c) Johna Smith, 1996
//
//  Given: z1,z2 - complex values
//  z1=a1+i*b1, z2=a2+i*b2
//  z1*z2=(a1*a2-b1*b2)+i*(a1*b2+a2*b1)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
struct complex
{
  float re;
  float im;
};
void show_complex(complex c) // this functions displays complex value
{
  printf("%f%+fi",c.re,c.im);
}
complex Mul(complex a, complex b)
{
  complex c;
  c.re=a.re*b.re-a.im*b.im;
  c.im=a.im*b.re+a.re*b.im;
  return c;
}
complex a={2,3},b={-1,2};
void main(void)
{
  show_complex(a);
  printf(" * ");
  show_complex(b);
  printf(" = ");
  show_complex(Mul(a,b));
}

Деление в экспоненциальной форме
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (division in M*exp(i*phi) form)
//  (c) Johna Smith, 1996
//
//  Given: z1,z2 - complex values
//  z1=M1*exp(i*phi1), z2=M2*exp(i*phi2)
//  z1/z2=M1/M2*exp(i*(phi1-phi2))
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi      3.1415926536
struct exp_complex
{
  float M;
  float phi;
};
void show_exp_complex(exp_complex c) // this functions displays complex value
{
  printf("%f*exp(%f*i)",c.M,c.phi);
}
exp_complex Div(exp_complex a,exp_complex b)
{
  exp_complex c;
  c.M=a.M/b.M;
  c.phi=a.phi-b.phi;
  return c;
}
exp_complex a={1,0}, b={-1,-Pi/2};
void main(void)
{
  show_exp_complex(a);
  printf(" / ");
  show_exp_complex(b);
  printf(" = ");
  show_exp_complex(Div(a,b));
}

Умножение в экспоненциальной форме
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
//////////////////////////////////////////////////////////////////////////////
//
//  Complex values operations (multiplication in M*exp(i*phi) form)
//  (c) Johna Smith, 1996
//
//  Given: z1,z2 - complex values
//  z1=M1*exp(i*phi1), z2=M2*exp(i*phi2)
//  z1*z2=M1*M2*exp(i*(phi1+phi2))
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define Pi      3.1415926536
struct exp_complex
{
  float M;
  float phi;
};
void show_exp_complex(exp_complex c) // this functions displays complex value
{
  printf("%f*exp(%f*i)",c.M,c.phi);
}
exp_complex Mul(exp_complex a,exp_complex b)
{
  exp_complex c;
  c.M=a.M*b.M;
  c.phi=a.phi+b.phi;
  return c;
}
exp_complex a={-1,Pi/2}, b={-1,-Pi/2};
void main(void)
{
  show_exp_complex(a);
  printf(" * ");
  show_exp_complex(b);
  printf(" = ");
  show_exp_complex(Mul(a,b));
}


Математическая статистика

Интервальная оценка среднего (рассмотрены случаи значений доверительной вероятности 0.95 и 0.99)
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating <x>, (x-dx,x+dx)
//  (c) Johna Smith, 1996
//
//  Method description:
//    dx=t(alpha,n)*Sn/sqrt(n)
//    <x>=Sum x/n
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N   6 // number of experiments
float x[N]={3,5,4,5,3,4}; // experimental data
double Student095(unsigned int n)
{
  double result;
  result=1.96+2.4/(n-1)+5.901610*pow(n-1,-2.372);
  return result;
}
double Student099(unsigned int n)
{
  double result;
  result=2.576+5.0/(n-1)+29.12178*pow(n-1,-2.591843);
  return result;
}
// this function calculates average x
double Average(void)
{
  double sum=0;
  for (int i=0;i<N;i++) sum+=x[i];
  return sum/N;
}
// this function calculates delta x
double DeltaX(void)
{
  double ave,sum=0,Sn;
  ave=Average();
  for (int i=0;i<N;i++) sum+=x[i]*x[i];
  Sn=sqrt(((double)N/(N-1))*((double)sum/N-ave*ave));
  return Student095(N)*Sn/sqrt(N); // you can use Student099 if you want
}
void main(void)
{
  printf("Data: ");
  for (int i=0;i<N;i++) printf("%f ",x[i]);
  printf("\nx=%f +/- %f",Average(),DeltaX());
}

Вычисление коэффициента Стьюдента
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating Student's coefficients t(alpha,n), if alpha=0.95 or 0.99
//  (c) Johna Smith, 1996
//
//  Method description:
//   To calculate Student's coefficients we use the following
//  approximation:
//    t(n)=C0+C1/(n-1)+C2(n-1)^(-C3)
//  Coefficients Ci for alpha=0.95 and alpha=0.99 is given
//  n - number of experiments
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double Student095(unsigned int n)
{
  double result;
  result=1.96+2.4/(n-1)+5.901610*pow(n-1,-2.372);
  return result;
}
double Student099(unsigned int n)
{
  double result;
  result=2.576+5.0/(n-1)+29.12178*pow(n-1,-2.591843);
  return result;
}
void main(void)
{
  printf("alpha=0.95 n=5 t5= %g\n",Student095(5));
  printf("alpha=0.99 n=15 t15= %g",Student099(15));
}

Теоретический коэффициент корреляции для двух величин
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating theoretical corellation coefficient for two values
//  (c) Johna Smith, 1996
//
//  Method description:
//    Corellation coefficient shows dependencies between two values
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N       3 // number of experiments
double x[N]={1,2,3};
double y[N]={2,4,5};
double TheoreticalCorellation(void)
{
  double sxy=0,sx2=0,sy2=0;
  for (int i=0;i<N;i++)
  {
    sxy+=x[i]*y[i];
    sx2+=x[i]*x[i];
    sy2+=y[i]*y[i];
  }
  return sxy/sqrt(sx2*sy2);
}
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%f %f\n",x[i],y[i]);
  printf("\nCorellation coefficient r= %f\n",TheoreticalCorellation());
}

Сглаживание экспериментальных данных по 3 точкам
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
//////////////////////////////////////////////////////////////////////////////
//
//  Smoothing experimental data by 3 points
//  (c) Johna Smith, 1996
//
//  Method description:
//    We smooth experimental data to enclose them to real curve.
//    All formulas use smallest sqares method to approximate data.
//    This program uses smoothing using 3-points groups
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N 10      // number of experiments
float data[N]={0.9,2.12,2.92,4.15,4.9,6.1,6.92,8.15,9.05,9.8};
float smoothed[N];
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%.4f ",data[i]);
  //smoothing
  smoothed[0]=(5*data[0]+2*data[1]-data[2])/6;
  for (i=1;i<N-1;i++) smoothed[i]=(data[i-1]+data[i]+data[i+1])/3;
  smoothed[N-1]=(5*data[N-1]+2*data[N-2]-data[N-3])/6;
  printf("\nSmoothed data:\n");
  for (i=0;i<N;i++) printf("%.4f ",smoothed[i]);
}

Сглаживание экспериментальных данных по 5 точкам
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
//////////////////////////////////////////////////////////////////////////////
//
//  Smoothing experimental data by 5 points
//  (c) Johna Smith, 1996
//
//  Method description:
//    We smooth experimental data to enclose them to real curve.
//    All formulas use smallest sqares method to approximate data.
//    This program uses smoothing using 5-points groups
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N  10      // number of experiments
float data[N]={0.9,2.12,2.92,4.15,4.9,6.1,6.92,8.15,9.05,9.8};
float smoothed[N];
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%.4f ",data[i]);
  //smoothing
  smoothed[0]=(3*data[0]+2*data[1]+data[2]-data[4])/5;
  smoothed[1]=(4*data[0]+3*data[1]+2*data[2]+data[3])/10;
  for (i=2;i<N-2;i++) smoothed[i]=(data[i-2]+data[i-1]+data[i]+
                                   data[i+1]+data[i+2])/5;
  smoothed[N-2]=(data[N-4]+2*data[N-3]+3*data[N-2]+4*data[N-1])/10;
  smoothed[N-1]=(3*data[N-1]+2*data[N-2]+data[N-3]-data[N-5])/5;
  printf("\nSmoothed data:\n");
  for (i=0;i<N;i++) printf("%.4f ",smoothed[i]);
}

Сглаживание экспериментальных данных по 7 точкам
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
//////////////////////////////////////////////////////////////////////////////
//
//  Smoothing experimental data by 7 points
//  (c) Johna Smith, 1996
//
//  Method description:
//    We smooth experimental data to enclose them to real curve.
//    All formulas use smallest sqares method to approximate data.
//    This program uses smoothing using 7-points groups
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N 10      // number of experiments
float data[N]={0.7,1.08,1.39,1.64,1.76,1.99,2.04,2.22,2.28,2.42};
float smoothed[N];
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%.4f ",data[i]);
  //smoothing
  smoothed[0]=(39*data[0]+8*data[1]-4*(data[2]+data[3]-data[5]) + data[4]-2*data[6])/42;
  smoothed[1]=(8*data[0]+19*data[1]+16*data[2]+6*data[3]-4*data[4] - 7*data[5]+4*data[7])/42;
  smoothed[2]=(-4*data[0]+16*data[1]+19*data[2]+12*data[3]+2*data[4] - 4*data[5]+data[6])/42;
  for (i=3;i<N-3;i++)
    smoothed[i]=(7*data[i]+6*(data[i-1]+data[i+1]) + 
          3*(data[i-2]+data[i+2])-2*(data[i-3]+data[i+3]))/21;
  smoothed[N-3]=(data[N-7]-4*data[N-6]+2*data[N-5]+12*data[N-4] + 
          9*data[N-3]+16*data[N-2]-4*data[N-1])/42;
  smoothed[N-2]=(4*data[N-7]-7*data[N-6]-4*data[N-5]+6*data[N-4] + 
          16*data[N-3]+19*data[N-2]+8*data[N-1])/42;
  smoothed[N-1]=(-2*data[N-7]+4*data[N-6]+data[N-5]-4*data[N-4] - 
          4*data[N-3]+8*data[N-2]+39*data[N-1])/42;
  printf("\nSmoothed data:\n");
  for (i=0;i<N;i++) printf("%.4f ",smoothed[i]);
}

Эмпирический коэффициент корреляции для двух величин
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating empirical corellation coefficient for two values
//  (c) Johna Smith, 1996
//
//  Method description:
//    Corellation coefficient shows dependencies between two values
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N 3 // number of experiments
double x[N]={1,2,3};
double y[N]={2,4,5};
double EmpiricalCorellation(void)
{
  double sxy=0,sx2=0,sy2=0,sx=0,sy=0;
  for (int i=0;i<N;i++)
  {
    sxy+=x[i]*y[i];
    sx2+=x[i]*x[i];
    sy2+=y[i]*y[i];
    sx+=x[i];
    sy+=y[i];
  }
  return (sxy-sx*sy/N)/(sqrt(sx2-sx*sx/N)*sqrt(sy2-sy*sy/N));
}
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%f %f\n",x[i],y[i]);
  printf("\nCorellation coefficient r= %f\n",EmpiricalCorellation());
}

Эмпирический коэффициент корреляции для трех величин
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating corellation coefficient for three values
//  (c) Johna Smith, 1996
//
//  Method description:
//    Corellation coefficient shows dependencies between three values
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N   3 // number of experiments
double x[N]={1,2,3};
double y[N]={2,4,5};
double z[N]={1,7,0};
double EmpiricalCorellation2D(double *i, double *j)
{
  double sij=0,si2=0,sj2=0,si=0,sj=0;
  for (int a=0;a<N;a++)
  {
    sij+=*(i+a)**(j+a);
    si2+=*(i+a)**(i+a);
    sj2+=*(j+a)**(j+a);
    si+=*(i+a);
    sj+=*(j+a);
  }
  return (sij-si*sj/N)/(sqrt(si2-si*si/N)*sqrt(sj2-sj*sj/N));
}
double Corellation3D(void)
{
  double rxz,ryz,rxy;
  // calculating partial corellation coefficients
  rxy=EmpiricalCorellation2D(x,y);
  rxz=EmpiricalCorellation2D(x,z);
  ryz=EmpiricalCorellation2D(y,z);
  return sqrt((rxz*rxz-ryz*ryz-2*rxy*rxz*ryz)/(1-rxy*rxy));
}
void main(void)
{
  printf("Data:\n");
  for (int i=0;i<N;i++) printf("%f %f %f\n",x[i],y[i],z[i]);
  printf("\nCorellation coefficient R= %f\n",Corellation3D());
}

Частный коэффициент корреляции для трех величин
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating partial corellation coefficient for three values
//  (c) Johna Smith, 1996
//
//  Method description:
//   Partial corellation coefficient shows dependencies between two values
//   using set of three values (we exclude influence of one of these values)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
// excluding influence of z, recieving corellation
// coefficients between x,y  y,z and x,z
double PartialCorellation3D(double rxy,double rxz,double ryz)
{
  return (rxy-rxz*ryz)/sqrt((1-rxz*rxz)*(1-ryz*ryz));
}
void main(void)
{
  double rxy=0.998753,rxz=0.69753,ryz=0.7473;
  printf("Data:\n");
  printf("%f %f %f\n",rxy,rxz,ryz);
  printf("\nPartial corellation coefficient r= %.10f\n",
  PartialCorellation3D(rxy,rxz,ryz));
}


Добавлено через 8 минут
Интегрирование функций

Метод трапеций
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
//////////////////////////////////////////////////////////////////////////////
//
//  Integrating function (trapecias method)
//  (c) Johna Smith, 1996
//
//  Method description:
//   We approximate area under function with trapecias, using
//   the same intervals by the X axe.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
double f(double x)
{
  // !!! you should specify your own function here
  return x;
}
double Integrate(double x1, double x2, double step)
{
  double x=x1;
  double I=0; // integral value
  while (x<x2-step)
  {
    // calculating trapecia's area as f(x)+f(x+step) * halfheight of it
    I += (f(x)+f(x+step))*step/2;
    x += step;
  }
  return I;
}
void main(void)
{
  printf("I=%f",Integrate(0,1,1e-4));
}

Метод Симпсона
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
//////////////////////////////////////////////////////////////////////////////
//
//  Integrating function (Simpson's method)
//  (c) Johna Smith, 1996
//
//  Method description:
//    We approximate area under function with parabolas, using
//    three points. We repeat calculation with less step if
//    the diffference between I(h) and I(h/2) is greater than eps
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
double f(double x)
{
  // you should specify your own function here
  return x*x*sin(x);
}
double Integrate(double x1, double x2, double step, double precision)
{
  double x=x1;
  double I1=0,I2=0; // integral value
  double s;
  do
  {
    I1=I2;
    s=(f(x1)-f(x2))/2;
    x=x1+step;
    while (x<x2)
    {
      s+=2*f(x)+f(x+step);
      x+=2*step;
    }
    I2=2*step*s/3;
    step/=2.0;  // try once more using less step
  } while (fabs(I1-I2)>precision);
  return I2;
}
void main(void)
{
  printf("I=%f",Integrate(0,1,1e-4,1e-6));
}

Метод прямоугольников
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
//////////////////////////////////////////////////////////////////////////////
//
//  Integrating function (rectanges method)
//  (c) Johna Smith, 1996
//
//  Method description:
//    We approximate area under function with rectangles, using
//    the same intervals by the X axe.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
double f(double x)
{
  // you should specify your own function here
  return x;
}
double Integrate(double x1, double x2, double step)
{
  double x=x1;
  double I=0; // integral value
  
  while (x<x2-step)
  {
    // calculating function at the halfpoint of interval
    I+=step*f(x+step/2);
    x+=step;
  }
  
  return I;
}
void main(void)
{
  printf("I=%f",Integrate(0,1,1e-4));
}

Метод Монте-Карло
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
//////////////////////////////////////////////////////////////////////////////
//
//  Integrating function (Monte-Karlo method)
//  (c) Johna Smith, 1996
//
//  Method description:
//   We generate random numbers from interval [a;b] where we want to
//   integrate our function and sum function values in this random points.
//   Here we use approximation of the integral:
//   I=<y>*sigma, where <y> is average value of the function f on [a;b]
//   and sigma is volume of integrating area (if we integrate by interval it
//   is its length)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
double f(double x)
{
  // you should specify your own function here
  return x*x*sin(x);
}
double Integrate(double x1, double x2, long int steps)
{
  double I=0; // integral value
  randomize();  // initializing random numbers generator
  for (long int i=0;i<steps;i++)
    I+=f((random(10000)/10000.0)*(x2-x1)+x1);
  I*=(x2-x1)/steps;
  return I;
}
void main(void)
{
  printf("I=%f",Integrate(0,1,1e5));
}


Календарь

Число дней между двумя датами
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating number of days between two dates
//  (c) Johna Smith, 1996
//
//  Method description:
//    We calculte two Julian dates and substract them
//    Result will be number of days between two dates
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
enum Month {January=1,February,March,April,May,June,July,August,
            September,October,November,December};
struct Date
{
  float day;
  Month month;
  int year;
};
double JulianDate(Date date)
{
  long int A,B,Y,M;
  if (date.month==January || date.month==February)
  {
    Y=date.year-1;
    M=date.month+12;
  } else
  {
    Y=date.year;
    M=date.month;
  }
  if ((date.year>1581) ||
      (date.year==1581 && date.month>October) ||
      (date.year==1581 && date.month==October && date.day>=15))
  {
    A=Y/100;
    B=2-A+A/4;
  } else B=0;
  return B+(long int)(365.25*Y)+(long int)(30.6001*(M+1))+date.day+1720994.5;
}
Date Revolution={17.5,October,1917};
Date Today={11,October,1996};
void main(void)
{
  printf("Days between %f/%d/%d and %f/%d/%d = %f",
         Revolution.day,Revolution.month,Revolution.year,
         Today.day,Today.month,Today.year,
         JulianDate(Today)-JulianDate(Revolution));
}

Определение дня недели по дате
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
//////////////////////////////////////////////////////////////////////////////
//
//  Calculating day of week by date
//  (c) Johna Smith, 1996
//
//  Method description:
//  Here we use special formulas different for January & February
//  and other monthes. Use this algorithm only for dates after
//  October, 15, 1528
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
enum Month {January=1,February,March,April,May,June,July,August,
            September,October,November,December};
struct Date
{
  float day;
  Month month;
  int year;
};
int DayOfWeek(Date date)
{
  float F;
  if (date.month<March)
    F=365*date.year+date.day+31*(date.month-1)+
      (int)((date.year-1)/4)-(int)(3*(int)((date.year-1)/100+1)/4);
  else
    F=365*date.year+date.day+31*(date.month-1)-(int)(0.4*date.month+2.3)+
      (int)(date.year/4)-(int)(3*(int)(date.year/100+1)/4);
  return  F-7*(int)(F/7)-1;
}
Date Today={11,October,1996};
void main(void)
{
  printf("%f/%d/%d is ",Today.day,Today.month,Today.year);
  switch (DayOfWeek(Today))
  {
    case -1:printf("Sunday\n");break;
    case 0:printf("Monday\n");break;
    case 1:printf("Tuesday\n");break;
    case 2:printf("Wednesday\n");break;
    case 3:printf("Thursday\n");break;
    case 4:printf("Friday\n");break;
    case 5:printf("Saturday\n");break;
  }
}

Вычисление даты по юлианскому календарю
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//////////////////////////////////////////////////////////////////////////////
//
//  Julian date calculating
//  (c) Johna Smith, 1996
//
//  Method description:
//    Julian date is the number of days since afternoon (GMT) of
//    January, 1, 4713 BC.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
enum Month {January=1,February,March,April,May,June,July,August,
            September,October,November,December};
struct Date
{
  float day;
  Month month;
  int year;
};
double JulianDate(Date date)
{
  long int A,B,Y,M;
  if (date.month==January || date.month==February)
  {
    Y=date.year-1;
    M=date.month+12;
  } else
  {
    Y=date.year;
    M=date.month;
  }
  if ((date.year>1581) ||
      (date.year==1581 && date.month>October) ||
      (date.year==1581 && date.month==October && date.day>=15))
  {
    A=Y/100;
    B=2-A+A/4;
  } 
  else 
    B=0;
  return B+(long int)(365.25*Y)+(long int)(30.6001*(M+1))+date.day+1720994.5;
}
Date Revolution={17.5,October,1917};
void main(void)
{
  printf("%f/%d/%d = %f",Revolution.day,Revolution.month,
                         Revolution.year,JulianDate(Revolution));
}


Решение систем линейных уравнений

Метод квадратного корня
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
//////////////////////////////////////////////////////////////////////////////
//
//  Solving system of linear equations (square root method)
//  (c) Johna Smith, 1996
//
//  Method description:
//    From given matrix A we build two auxulary triangular matrixes S & D
//    Then solving equations Gy=b (where g[I][j]=s[j][I]*d[j]) and Sx=y
//    we obtain vector x.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N    4     // size of matrix
#define N1   N+1
float matrix[N][N1]=
    {{10.9,1.2,2.1,0.9,23.2},
  {1.2,11.2,1.5,2.5,38.1},
  {2.7,1.5,9.8,1.3,40.3},
  {0.9,2.5,1.3,12.1,58.2}
    };
void ShowMatrix(void)
{
  for (int i=0;i<N;i++)
  {
    for (int j=0;j<N;j++)
      printf("%+f*x%d",matrix[i][j],i+1);
    printf("=%f\n",matrix[i][N]);
  }
}
void main(void)
{
  // Variables declaration
  register char i,j,k;
  float s[N][N],x[N],y[N],d[N],sum;
  // Printing given matrix
  ShowMatrix();
  // Building matrixes S and D
  for (j=0;j<N;j++)
    for (i=0;i<=j;i++)
    {
      sum=matrix[i][j];
      for (k=0;k<i;k++)
        sum-=s[k][i]*d[k]*s[k][j];
      if (i==j)
      {
        d[i]=(sum>0?1:-1);
        s[i][j]=sqrt(fabs(sum));
      }
      else s[i][j]=sum/(d[i]*s[i][i]);
    }
  // Solving equation Gy=b (G: g[I][j]=s[j][I]*d[j])
  for (i=0;i<N;i++)
  {
    y[i]=matrix[i][N];
    for (j=0;j<i;j++) y[i]-=s[j][i]*d[j]*y[j];
    y[i]/=s[i][i]*d[i];
  }
  // Solving equation Sx=y
  for (i=N-1;i>=0;i--)
  {
    x[i]=y[i];
    for (j=i+1;j<N;j++) x[i]-=s[i][j]*x[j];
    x[i]/=s[i][i];
  }
  // Printing solution
  printf("\nSolution:\n");
  for (i=0;i<N;i++)
  printf("x%d=%f\n",i+1,x[i]);
}

Метод Некрасова
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
//////////////////////////////////////////////////////////////////////////////
//
//  Solving system of linear equations (Nekrasov method)
//  (c) Johna Smith, 1996
//
//  Method description:
//    This is an iterative method of solving systems of linear equations.
//    Criteria for end of iterations is ||xn||-||x(n-1)||<epsilon, where
//    ||x|| is norm of vector x.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N     4     // size of matrix
#define N1    N+1
float matrix[N][N1]=
      {{10.9,1.2,2.1,0.9,23.2},
       {1.2,11.2,1.5,2.5,38.1},
       {2.7,1.5,9.8,1.3,40.3},
       {0.9,2.5,1.3,12.1,58.2}
     };
void ShowMatrix(void)
{
  for (int i=0;i<N;i++)
  {
    for (int j=0;j<N;j++)
      printf("%+f*x%d",matrix[i][j],i+1);
    printf("=%f\n",matrix[i][N]);
  }
}
// this function calculates ||x||
float norma(float *x)
{
  float sum=0;
  for (unsigned i=0;i<N;i++) sum+=fabs(*(x+i));
  return sum;
}
void main(void)
{
  // Variables declaration
  register char i,j;
  float x_current[N],x_next[N],sum1,sum2,epsilon;
  epsilon=1e-7;
  // Printing given matrix
  ShowMatrix();
  // Solving
  do
  {
    for (i=0;i<N;i++) x_current[i]=x_next[i];
    for (i=0;i<N;i++)
    {
      sum1=0;sum2=0;
      for(j=0;j<i;j++) sum1+=matrix[i][j]*x_next[j];
      for(j=i+1;j<N;j++) sum2+=matrix[i][j]*x_current[j];
      x_next[i]=(matrix[i][N]-sum1-sum2)/matrix[i][i];
    }
  } while (fabs(norma(x_current)-norma(x_next))>epsilon);
  // Printing solution
  printf("\nSolution:\n");
  for (i=0;i<N;i++)
    printf("x%d=%f\n",i+1,x_current[i]);
}

Метод 3-диагональных матриц
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
//////////////////////////////////////////////////////////////////////////////
//
//  Solving system of linear equations (3-diagonal matrix)
//  (c) Johna Smith, 1996
//
//  Method description:
//    This method is developed for 3-diagonal matrixes:
//    (a11 a12 0 ..............)
//    (a21 a22 a23 0 ..........)
//    (0   a32 a33 a34 0 ......) = M
//    (........................)
//    (.......... 0 an(n-1) ann)
//    Cause all elements except three diagonals is 0 we'll store
//    only non zero elements in the A,B,C arrays. Array D stores
//    vector B (M*X=B). We build two auxulary vectors P and Q and
//    then calculate vector X using the following recurrent formula:
//      X(k-1)=P(k)*X(k)+Q(k)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#define N    20
#define N1   N+1
void main(void)
{
  // Variables declaration
  register int i,i1,k;  // using processor registers for better performance
  float a[N]={0,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1};
  float b[N]={2,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,2};
  float c[N]={-1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0};
  float d[N]={0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0};
  float x[N1],p[N1],q[N1],tmp;
  // Some initializations
  p[0]=0;
  q[0]=0;
  x[N]=0;
  c[N-1]=0; // to avoid overflow
  // Building vectors P and Q
  for (i=0;i<N;i++)
  {
    i1=i+1;
    if ((tmp=b[i]+a[i]*p[i])==0)
    {
      printf("Wrong method for such system...\n");
      return;
    };
    p[i1]=-c[i]/tmp;
    q[i1]=(d[i]-a[i]*q[i])/tmp;
  }
  // Building vector X
  for (k=N;k>0;k--) x[k-1]=p[k]*x[k]+q[k];
  // Printing solution
  printf("Solution:\n");
  for (i=0;i<N;i++)
    printf("x%d=%f\n",i+1,x[i]);
}

Метод Гаусса
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
//////////////////////////////////////////////////////////////////////////////
//
//  Solving system of linear equations (Gauss method)
//  (c) Johna Smith, 1996
//
//  Method description:
//    This methos based on excluding variables from equations.
//    Using first equation we exclude x1 from all other equations
//    Using second equation equation we exclude x2 from all equations
//   (excluding first) and so on. So we have triangular matrix from which
//    we can easily get vector X. (AX=B)
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N     4     // size of matrix
#define N1    N+1
float matrix[N][N1]=
                  {{2,5,4,1,28},
                   {1,3,2,1,17},
                   {2,10,9,7,77},
                   {3,8,9,2,54}
                  };
void ShowMatrix(void)
{
  for (int i=0;i<N;i++)
  {
    for (int j=0;j<N;j++)
      printf("%+f*x%d",matrix[i][j],i+1);
    printf("=%f\n",matrix[i][N]);
  }
}
void main(void)
{
  // Variables declaration
  float tmp,x[N1];
  register short int i,j,k;
  // Printing given matrix
  ShowMatrix();
  // Main loop
  for (i=0;i<N;i++)
  {
    // Excluding variable x[i] from equations
    tmp=matrix[i][i];
    for (j=N;j>=i;j--) matrix[i][j]/=tmp;
    for (j=i+1;j<N;j++)
    {
      tmp=matrix[j][i];
      for (k=N;k>=i;k--)
        matrix[j][k]-=tmp*matrix[i][k];
    }
  }
  // Calculating vector x
  x[N-1]=matrix[N-1][N];
  for (i=N-2;i>=0;i--)
  {
    x[i]=matrix[i][N];
    for (j=i+1;j<N;j++) x[i]-=matrix[i][j]*x[j];
  }
  // Printing solution
  printf("\nSolution:\n");
  for (i=0;i<N;i++)
    printf("x%d=%f\n",i+1,x[i]);
}

Метод Гаусса с выбором максимального элемента
Кликните здесь для просмотра всего текста
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
//////////////////////////////////////////////////////////////////////////////
//
//  Solving system of linear equations (Gauss method with max element selecting)
//  (c) Johna Smith, 1996
//
//  Method description:
//    This methos based on excluding variables from equations.
//    Using first equation we exclude x1 from all other equations
//    Using second equation equation we exclude x2 from all equations
//    (excluding first) and so on. So we have triangular matrix from which
//    we can easily get vector X. (AX=B)
//    'Selecting max element' means that on each step we select equation
//    with max element in current column and use it as next equation to exclude
//    variable.
//
//////////////////////////////////////////////////////////////////////////////
#include <stdio.h>
#include <math.h>
#define N      4     // size of matrix
#define N1     N+1
float matrix[N][N1]=
    {{2,5,4,1,28},
  {1,3,2,1,17},
   {2,10,9,7,77},
  {3,8,9,2,54}
        };
void ShowMatrix(void)
{
  for (int i=0;i<N;i++)
  {
    for (int j=0;j<N;j++)
      printf("%+f*x%d",matrix[i][j],i+1);
    printf("=%f\n",matrix[i][N]);
  }
}
void main(void)
{
  // Variables declaration
  float max,swp,tmp,x[N1];
  register short int row_with_max_element,i,j,k;
  // Printing given matrix
  ShowMatrix();
  // Main loop
  for (i=0;i<N;i++)
  {
    // Searching for max element in the current column (i)
    max=fabs(matrix[i][i]);
    row_with_max_element=i;
    for (j=i+1;j<N;j++)
      if (fabs(matrix[j][i])>max)
      {
        row_with_max_element=j;
        max=fabs(matrix[j][i]);
      }
    // Swapping 2 lines of matrix - row_with_max_element & i
    if (row_with_max_element!=i)
    for (j=0;j<N1;j++)
    {
      swp=matrix[row_with_max_element][j];
      matrix[row_with_max_element][j]=matrix[i][j];
      matrix[i][j]=swp;
    }
    // Excluding variable x[i] from equations
    tmp=matrix[i][i];
    for (j=N;j>=i;j--) matrix[i][j]/=tmp;
    for (j=