Коллекция алгоритмов, в помощь начинающим |
Здравствуйте, гость ( Авторизация | Регистрация )
Коллекция алгоритмов, в помощь начинающим |
4.01.2010 - 14:21
Сообщение
#1
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Выложенные здесь алгоритмы преследуют исключительно учебные цели.
Код неоптимизирован, местами морально устарел и показывает только принцип решения той или иной задачи. Однако он вполне рабочий и может быть использован с соответствующей дорботкой под Ваши собственные нужды. Буду рада, если не только я буду их пополнять) Содержание: Перевод из ... в ... Сортировка массивов Графика Матрицы Специальные функции Методы поиска Решение систем линейных уравнений Решение нелинейных уравнений Решение диф. уравнений Будет еще дополнено. -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
|
4.01.2010 - 14:29
Сообщение
#2
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Градусы, минуты, секунды - в градусы:
» Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Градусы - в градусы, минуты, секунды: » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Комплексных величин - в экспоненциальную форму » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Градусы по Фаренгейту - в градусы по Цельсию » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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); } Полярные - в декартовы » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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); } Сферические - в декартовы » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
4.01.2010 - 14:37
Сообщение
#3
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Сортировка массивов
Binary insertions » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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(); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
5.01.2010 - 13:51
Сообщение
#4
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Графика
Фрактал (листья папоротника) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Рисование линии (по Брезенхэму) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Рисование окружности (по Брезенхэму) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Сглаживание кривой В-сплайном » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Рисование куба » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Рисование тора » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } Рисование 3D объекта » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 объекта » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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; } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
10.01.2010 - 19:21
Сообщение
#5
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Матрицы:
Определитель » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); Обращение » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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]); } Умножение » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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); } Транспонирование » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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]); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
12.01.2010 - 13:33
Сообщение
#6
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Специальные функции:
arccos x (x - real) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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) » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } arcsin x » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Функция Бесселя » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Гамма-функция » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } Бета-функция » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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)); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
14.01.2010 - 18:16
Сообщение
#7
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Методы поиска
Binary search in sorted string » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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 » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
22.01.2010 - 17:05
Сообщение
#8
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Решение систем линейных уравнений
Метод квадратного корня » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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]); } Метод Некрасова » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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-диагональных матриц » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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]); } Метод Гаусса » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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]); } Метод Гаусса с выбором максимального элемента » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // 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=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]); } Метод итераций » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving system of linear equations (iterations method) // (c) Johna Smith, 1996 // // Method description: // This program automatically transforms system of equations to // iterative form and solves it. // Iterative form: // x1=a11*x1+a12*x2+... + a1n*xn+b1 // x2=a21*x1+a22*x2+... + a2n*xn+b2 // ... // xn=an1*x1+an2*x2+... + ann*xn+bn // Setting first iteration of vector X we can get next iteration from // such form of system of equations. Iterations will be finished when // difference between two consequitive iterations will be less than epsilon // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> #define N 3 // size of matrix #define N1 N+1 float matrix[N][N1]= {{14.38,-2.41,1.39,5.86}, {1.84,25.36,-3.31,-2.28}, {2.46,-3.49,16.37,4.47} }; int maxiterations=10; // maximum number of iterations float epsilon=0.0001; // required accuracy 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 x[N],y[N],t; register short int i,j,k; int iterations=0; // Printing given matrix ShowMatrix(); // setting first iteration of vector X for (i=0;i<N;i++) x[i]=matrix[i][N]; do { for (i=0;i<N;i++) { t=-matrix[i][N]; for (j=0;j<N;j++) t+=matrix[i][j]*x[j]; y[i]=(-t+matrix[i][i]*x[i])/matrix[i][i]; } iterations++; k=0; // checking solution while (fabs(x[k]-y[k])<epsilon && k<N) k++; // new iteration becomes old for(i=0;i<N;i++) x[i]=y[i]; } while (k!=N && iterations<maxiterations); if (iterations==maxiterations) { printf("Iterations are very slow..."); } else { // Printing solution printf("\nSolution:\n"); for (i=0;i<N;i++) printf("x%d=%f\n",i+1,x[i]); printf("%d iterations were made",iterations); } } Метод Монте-Карло » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving system of linear equations (Monte-Karlo method) // (c) Johna Smith, 1996 // // Method description: // We should translate given system of linear equations to special form: // x1=a(11)x1+a(12)x2+...+a(1n)xn+a(1n+1) // x2=a(21)x2+a(22)x2+...+a(2n)xn+a(2n+1) // .... // xn=a(n1)x1+a(n2)x2+...+a(nn)xn+a(nn+1) // n // where SUM |aij|<1 (i=1, 2, ... ,n) // j=1 // Lets divide interval [0;1] into N+1 smaller intervals. // Imagine a particle that is moving randomly at [0;1] interval. // Remember its moves until it reachs one of the bounds of the interval. // So we'll get a trajectory of the particle s(i),s(j),s(k),s(l),s(m),... // If particle moves from interval s(i) to s(j) we'll write v(ij) // y[i]=v(ij)*v(jk)*...*v(tm)*w(m), v(ij)=sign(aij), v(jk)=sign(ajk) // n // w(m)=a(mn+1)/(1-SUM a(mj)), x[i]=y/M, where M is number of particle runs // j=1 // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> #include <stdlib.h> #define N 3 // size of matrix #define N1 N+1 #define M 1000 // number of tries float matrix[N][N1]= {{3,1,-1,-2}, {1,4,-2,-3}, {2,-1,5,-15} }; void ShowMatrix(void) { for (int i=0;i<N;i++) { printf("x%d=",i+1); for (int j=0;j<N;j++) printf("%+f*x%d",matrix[i][j],j+1); printf("%+f\n",matrix[i][N]); } } void main(void) { // Variables declaration float b[N][N],w[N],y,c,x[N],l; register short int i,j,v; unsigned char finish; int row, // currently analysed row tries; // number of tries to calculate true root by generating randoms // Transforming matrix for (i=0;i<N;i++) { l=0; for (j=0;j<N1;j++) l+=fabs(matrix[i][j]); v=(matrix[i][i]>0?-1:1); for (j=0;j<N1;j++) matrix[i][j]=v*(matrix[i][j])/l+(i==j?1:0); } // Printing transformed matrix ShowMatrix(); // filling matrix B for (i=0;i<N;i++) { b[i][0]=fabs(matrix[i][0]); for(j=1;j<N;j++) b[i][j]=b[i][j-1]+fabs(matrix[i][j]); } // filling matrix w for(i=0;i<N;i++) w[i]=matrix[i][N]/(1-b[i][N-1]); // for each equation in the system for(i=0;i<N;i++) { tries=0; row=i; y=0; v=1; while (tries<M) { // generating random number c=random(32767)/32767.0; j=N-1; finish=0; while(j>=0 && !finish) { if (c>b[row][j]) { if (j==N-1) { tries++; y+=v*w[row]; row=i;v=1; } else { v*=(matrix[row][j+1]>=0?1:-1); row=j+1; } finish=1; } j--; } if (!finish) { v*=(matrix[row][j+1]>=0?1:-1); row=0; } } // calculating root x[i]=y/M; } // Printing solution printf("\nSolution:\n"); for (i=0;i<N;i++) printf("x%d=%f\n",i+1,x[i]); } Метод ортогонализации » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving system of linear equations (orthogonalization method) // (c) Johna Smith, 1996 // // Method description: // Given system of linear equations: // a11 x1 + a12 x2 + ... + a1n xn + a(1n+1) = 0 // a21 x1 + a22 x2 + ... + a2n xn + a(2n+1) = 0 // // an1 x1 + an2 x2 + ... + ann xn + a(nn+1) = 0 // // Left part of each equation is a result of scalar multiplication of // two vectors: ai=(ai1,ai2,...,ain,a(in+1)) and x=(x1,x2,..,xn,1) // So to solve this system we need only to build vector x which will be // orthogonal to each of ai vectors. // n+1 // Let u1=a1, b1=u1/sqrt(SUM u(1j)^2) // j=1 // // Other ui we can get from the following iterative formula: // k // u(k+1)=u(k+1)-SUM {u(k+1), bj}bj, where {..} means scalar multiplication // j=1 // n+1 // b(k+1)=u(k+1)/sqrt(SUM u(k+1,j)^2) // j=1 // And finally we can obtain roots of this system from this formula: // // xi=b(n+1,i)/b(n+1n+1) // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> #include <stdlib.h> #define N 3 // size of matrix #define N1 N+1 double matrix[N1][N1]= {{3,1,-1,-2}, {1,4,-2,-3}, {2,-1,5,-15} }; 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=0\n",matrix[i][N]); } } void main(void) { // Variables declaration double c[N],x[N],s; int i,j,k,l,m; ShowMatrix(); // expanding system by vector (0,0,0,...,0,1) for (i=0;i<N;i++) matrix[N][i]=0; matrix[N][N]=1; // for all equations for (i=0;i<N1;i++) { if (i>0) { k=0; // make some iterations to increase accuracy of calculations while (k<=3) { for (l=0;l<i;l++) { c[l]=0; for(m=0;m<N1;m++) c[l]+=matrix[l][m]*matrix[i][m]; } for (j=0;j<N1;j++) { s=0; for(l=0;l<i;l++) s+=c[l]*matrix[l][j]; matrix[i][j]-=s; } k++; } } // normalizing vector s=0; for (k=0;k<N1;k++) s+=matrix[i][k]*matrix[i][k]; s=sqrt(s); for (j=0;j<N1;j++) matrix[i][j]/=s; } s=matrix[N][N]; for (j=0;j<N;j++) x[j]=matrix[N][j]/s; // Printing solution printf("\nSolution:\n"); for (i=0;i<N;i++) printf("x%d=%f\n",i+1,x[i]); } Метод Зэйделя » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving system of linear equations (Zaidel method) // (c) Johna Smith, 1996 // // Method description: // Given: system of linear equations // a11 x1 + a12 x2 + ... + a1n xn = b1 // a21 x1 + a22 x2 + ... + a2n xn = b2 // ... // an1 x1 + an2 x2 + ... + ann xn = bn // // We use following iterative formula to solve this system: // // 1 i-1 N // xi(j+1)=xi(j)- --- ( SUM aik xk(j+1) + SUM aik xk(j) -bi) // aii k=1 k=i // // where xi(j+1) is j+1-th iteration, xi(j) is j-th iteration // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> #define N 3 // size of matrix #define N1 N+1 float matrix[N][N1]= {{10,1,1,12}, {2,10,1,13}, {2,2,10,14} }; float z[N]={0,0,0}; // first iteration float epsilon=0.0001; // required accuracy 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 x[N]; short int i,j; int iterations=0,finish=0; // Printing given matrix ShowMatrix(); while (!finish) { finish=1; for (i=0;i<N;i++) { x[i]=-matrix[i][N]; for (j=0;j<N;j++) x[i]+=matrix[i][j]*z[j]; // don't stop iterations until required accuracy is reached if (fabs(x[i]/matrix[i][i])>=epsilon) finish=0; x[i]=z[i]-x[i]/matrix[i][i]; // next iteration z[i]=x[i]; } iterations++; } // Printing solution printf("\nSolution:\n"); for (i=0;i<N;i++) printf("x%d=%f\n",i+1,x[i]); printf("%d iterations were made",iterations); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
28.01.2010 - 19:18
Сообщение
#9
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Решение нелинейных уравнений
Метод хорд » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (chords method) // (c) Johna Smith, 1996 // // Method description: // given: f(x)=0, [a;b], f(a)*f(b)<0 // find x0: f(x0)=0 // We approximate curve by chord between its borders and get // point where the chord crosses X axe as a new border of // interval with the root. When interval becomes small enough // we can take one of its borders as a solution. // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> // this function returns value of f(x) double f(double x) { // sin 2x - ln x = 0 return sin(2*x)-log(x); } double Solve(double a, double b, double epsilon) { double u,v,x,y; u=f(a); v=f(b); do { // x - new border x=(a*v-b*u)/(v-u); y=f(x); // determine whether x is left or right border // f(a)*f(b)<0 condition must be true if (y*u<0) { b=x; v=y; } else { a=x; u=y; } } while (b-a>epsilon); return x; } void main(void) { printf("sin 2x - ln x = 0, [a;b]=[1.3;1.5]"); printf("\nx0=%f",Solve(1.3,1.5,1e-9)); } Метод дихотомии » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (dichotomia method) // (c) Johna Smith, 1996 // // Method description: // given: f(x)=0, [a;b], f(a)*f(b)<0 // find x0: f(x0)=0 // We split interval into two parts and select that part that // contains root of the equation by checking condition // f(a)*f(b)<0 We take middle point of the interval as a new border of // interval with the root. When interval becomes small enough // we can take one of its borders as a solution. // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> // this function returns value of f(x) double f(double x) { // sin 2x - ln x = 0 return sin(2*x)-log(x); } double Solve(double a, double b, double epsilon) { double x; do { // x - new border x=(b+a)/2; // determine whether x is left or right border // f(a)*f(b)<0 condition must be true if (f(x)*f(a)<0) b=x; else a=x; } while (b-a>epsilon); return x; } void main(void) { printf("sin 2x - ln x = 0, [a;b]=[1.3;1.5]"); printf("\nx0=%f",Solve(1.3,1.5,1e-9)); Метод Эйткена-Стеффенсона » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (Eitken-Steffenson method) // (c) Johna Smith, 1996 // // Method description: // This method is used for solving equations like x=f(x) // (see also iterations method) // Here we use the following iterations: // x(n+1)=(x0*x2-x1^2)/(x0-2x1+x2) // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> // this function returns value of f(x) double f(double x) { // x = x + 0.37(sin 2x - ln x) return x+0.37*(sin(2*x)-log(x)); } double Solve(double x, double epsilon) { double x1,x2,xold,tmp; do { xold=x; x1=f(x); x2=f(x1); tmp=x-2*x1+x2; x=(x*x2-x1*x1)/tmp; } while (fabs(xold-x)>epsilon && tmp!=0); return x; } void main(void) { printf("sin 2x - ln x = 0, [a;b]=[1.3;1.5]"); printf("\nx0=%f",Solve(1.4,1e-6)); } Метод итераций » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (iterations method) // (c) Johna Smith, 1996 // // Method description: // This method is used for solving equations like x=f(x) // We should set first approximation of the root and then // make some iterations x=f(x). When the difference between // two consequental iterations will be less than epsilon we // can say that we found the root. // To apply this method to equations like F(x)=0 we should // transform it into iterational form x=f(x), where f(x): // 1) defined on [a;b] // 2) for every point on [a;b] exists f'(x) // 3) for every x from [a;b] f(x) is in [a;b] // 4) exists number q : |f'(x)|<=q<1 for all x from [a;b] // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> // this function returns value of f(x) double f(double x) { // x = x + 0.37(sin 2x - ln x) return x+0.37*(sin(2*x)-log(x)); } double Solve(double x, double epsilon) { double x1; do { x1=x; x=f(x); } while (fabs(x1-x)>epsilon); return x; } void main(void) { printf("sin 2x - ln x = 0"); printf("\nx0=%f",Solve(1.4,1e-6)); } Метод Ньютона » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (Newton method) // (c) Johna Smith, 1996 // // Method description: // Here we use the following iterations: // x(n+1)=xn-f(xn)/f'(xn) // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> // this function returns value of f(x) double f(double x) { // sin 2x - ln x = 0 return sin(2*x)-log(x); } double Solve(double x, double epsilon) { double x1; do { x1=x; x=x-f(x)*epsilon/(f(x+epsilon)-f(x)); } while (fabs(x1-x)>epsilon); return x; } void main(void) { printf("sin 2x - ln x = 0, [a;b]=[1.3;1.5]"); printf("\nx0=%f",Solve(1.4,1e-9)); } Метод Лобачевского » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // Solving nonlinear equations (Lobachevsky method) // (c) Johna Smith, 1996 // // Method description: // Given: a0+a1x+a2x^2+...+anx^n=0 // This method allows to find modulus of the greatest root of this equation // even if it's complex. But in last case there can appear several messages // about impossibilty of calculation root of negative number. // The main idea of this method is to change given equation to other // equation which roots equals to powered roots of given equation. For example // if roots of the given equation are x0,x1,.. xn then roots of new equation // will be x0^2, x1^2, ..., xn^2. Repeating this operation we get an equation // where one root is much greater than other ones. So we can easily // obtain modulus of the greatrest root of the given equation. // To obtain other roots of equation we need to divide given equation // by (x-x0) (where x0 is found root) and apply this method to result. // ////////////////////////////////////////////////////////////////////////////// #include <stdio.h> #include <math.h> #define N 4 #define N1 N+1 #define Iterations 15 // number of iterations double a[N1]={24,-50,35,-10,1}; void main(void) { double r,b[N1],c[N1],g,bi,d; int z,k; // printing given equation printf("%f",a[0]); for(int i=1;i<N1;i++) printf("%+fx^%d",a[i],i); printf("=0\n\n"); // preparing auxiliary arrays b and c for (i=0;i<N1;i++) { b[i]=a[i]/a[N]; c[i]=0; } // setting required parameters r=1/2.0; g=1; // make all iterations for(int y=0;y<Iterations;y++) { // calculate coefficients c[i] (coefficients of new equation) z=1; for(i=0;i<N1;i++) { bi=z*b[i]; k=(i+1)/2; for(int j=i%2;j<N1;j+=2) { c[k]+=bi*b[j]; k++; } z=-z; } d=z*c[N-1]; // check whether we could calculate root of d if(d>0) { // calculating and printing new iteration g*=powl(d,r); printf("%f\n",g); for (i=0;i<N1;i++) { // preparing data for next iteration b[i]=c[i]/powl(d,N-i); c[i]=0; } b[N-1]=z; b[N]=-z; } else { // d is negative - can't calculate root for(i=0;i<N1;i++) { // preparing data for next iteration b[i]=c[i]; c[i]=0; } printf("no iteration (can't calculate root from negative number)\n"); } r/=2.0; } } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
|
5.03.2010 - 15:26
Сообщение
#10
|
|
Старейшина Текущее настроение: Вст. ник | Цитата Группа: Легенда Сообщений: 12414 Регистрация: 8.06.2009 Пользователь №: 27737 Из: сострадания к ближнему Награды: 205 Подарки: 171 Пол: ? Репутация: 2878 |
Решение дифференциальных уравнений
Метод Эйлера » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving differential equations (Eiler method) // (c) Johna Smith, 1996 // // Method description: // Given differential equation y'=f(x,y) and starting conditions // x=x0, y(x0)=y0. We should find solution of this equation // at [a;b] interval. // y(i+1) calculated as y(i) + delta y(i) // delta y=h*f(x,y(i)) // // In this example y'=cos(y)+3x, [a;b]=[0;1], x0=0, y0=1.3 // ////////////////////////////////////////////////////////////////////////////// #include <math.h> #include <stdio.h> const float a=0,b=1; // bound of the interval const int num_points=10; // number of point to solve float x0=0,y0=1.3; // initial conditions int M=1; float f(float x, float y) { return (cos(y)+3*x); // y'=cos(y)+3*x } void calculate(int m,float *y) { float x,yi,h; h=(b-a)/(num_points*m); yi=y0; x=x0; for (int i=0;i<num_points;i++) { for (int k=0;k<m;k++) { yi+=h*f(x,yi); x+=h; } *(y+i)=yi; } } void main(void) { float yh[num_points],yh2[num_points]; calculate(M,yh); calculate(2*M,yh2); // doubled step for better accuracy // epsilon is difference between solutions with single // and double steps printf("X\t\tYH\t\tYH2\t\tEPSILON\n"); for (int i=0;i<num_points;i++) printf("%f\t%f\t%f\t%f\n",(x0+((i+1)*(b-a))/num_points), yh[i],yh2[i],yh2[i]-yh[i]); } Метод Адамса » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving differential equations (Adams method) // (c) Johna Smith, 1996 // // Method description: // Given differential equation y'=f(x,y) and starting conditions // x=x0, y(x0)=y0. We should find solution of this equation // at [a;b] interval. // We use Runge-Kutta method to process first 'num_starting_points' // points and then use Adams extrapolation to process other points. // In this example y'=x+y, [a;b]=[0;2], x0=0, y0=1 // ////////////////////////////////////////////////////////////////////////////// #include <math.h> #include <stdio.h> const float a=0,b=2; // bounds of the interval const int num_points=10, // number of points to solve num_starting_points=4; // number of points to solve with Runge-Kutta method float x0=0,y0=1; // starting conditions float f(float x, float y) { return x+y; // y'=x+y } // this function realises Runge-Kutta method for n starting points void calculate(float *y) { float k1,k2,k3,k4,x,yi,h; h=(b-a)/num_points; // step yi=y0; x=x0; for (int i=0;i<num_starting_points;i++) { k1=h*f(x,yi); k2=h*f(x+h/2,yi+k1/2); k3=h*f(x+h/2,yi+k2/2); k4=h*f(x+h,yi+k3); yi+=(k1+2*k2+2*k3+k4)/6; x+=h; *(y+i+1)=yi; } } void main(void) { float y[num_points+1],h; y[0]=y0; // apply Runge-Kutta method calculate(y); h=(b-a)/num_points; // extrapolating for (int i=num_starting_points;i<num_points;i++) y[i] = y[i-1]+h/24*(55*f(x0+(i-1)*h,y[i-1])- 59*f(x0+(i-2)*h,y[i-2])+ 37*f(x0+(i-3)*h,y[i-3])- 9*f(x0+(i-4)*h,y[i-4])); printf("X\t\tY\t\tExact solution\n"); for (i=0;i<num_points;i++) printf("%f\t%f\t%f\n",(x0+i*h),y[i],(2*exp(x0+i*h)-(x0+i*h)-1)); } Метод Эйлера-Коши » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving differential equations (Eiler-Coshi method) // (c) Johna Smith, 1996 // // Method description: // Given differential equation y'=f(x,y) and starting conditions // x=x0, y(x0)=y0. We should find solution of this equation // at [a;b] interval. // y(i+1) calculated as y(i) + delta y(i) // delta y=h*(f(x,y(i))+f(x,z))/2, z=y(i)+h*(f(x(i),y(i)) // // In this example y'=cos(y)+3x, [a;b]=[0;1], x0=0, y0=1.3 // ////////////////////////////////////////////////////////////////////////////// #include <math.h> #include <stdio.h> const float a=0,b=1; // bound of the interval const int num_points=10; // number of point to solve float x0=0,y0=1.3; // initial conditions int M=1; float f(float x, float y) { return (cos(y)+3*x); // y'=cos(y)+3*x } void calculate(int m,float *y) { float x,yi,h,z; h=(b-a)/(num_points*m); yi=y0; x=x0; for (int i=0;i<num_points;i++) { for (int k=0;k<m;k++) { z=yi+h*f(x,yi); yi+=h*(f(x,yi)+f(x,z))/2; x+=h; } *(y+i)=yi; } } void main(void) { float yh[num_points],yh2[num_points]; calculate(M,yh); calculate(2*M,yh2); // doubled step for better accuracy // epsilon is difference between solutions with single // and double steps printf("X\t\tYH\t\tYH2\t\tEPSILON\n"); for (int i=0;i<num_points;i++) printf("%f\t%f\t%f\t%f\n",(x0+((i+1)*(b-a))/num_points), yh[i],yh2[i],yh2[i]-yh[i]); } Метод Рунге-Кутта » Кликните сюда для просмотра оффтоп текста.. « Код ////////////////////////////////////////////////////////////////////////////// // // Solving differential equations (Runge-Kutta method) // (c) Johna Smith, 1996 // // Method description: // Given differential equation y'=f(x,y) and starting conditions // x=x0, y(x0)=y0. We should find solution of this equation // at [a;b] interval. // y(i+1) calculated as y(i) + delta y(i) // delta y =1/6 k1 + 1/3 k2 + 1/3 k3 + 1/6 k4 // k1=hf(x,y) // k2=hf(x+h/2,y+k1/2) // k3=hf(x+h/2,y+k2/2) // k4=hf(x+h,y+k3) // // In this example y'=cos(x+y)+x-y, [a;b]=[0;1], x0=0, y0=0 // ////////////////////////////////////////////////////////////////////////////// #include <math.h> #include <stdio.h> const float a=0,b=1; // bound of the interval const int num_points=10; // number of point to solve float x0=0,y0=0; // initial conditions int M=1; float f(float x, float y) { return (cos(x+y)+x-y); // y'=cos(x+y)+x-y } void calculate(int m,float *y) { float k1,k2,k3,k4,x,yi,h; h=(b-a)/(num_points*m); yi=y0; x=x0; for (int i=0;i<num_points;i++) { for (int k=0;k<m;k++) { // calculatin coefficients k k1=h*f(x,yi); k2=h*f(x+h/2,yi+k1/2); k3=h*f(x+h/2,yi+k2/2); k4=h*f(x+h,yi+k3); yi+=(k1+2*k2+2*k3+k4)/6; x+=h; } *(y+i)=yi; } } void main(void) { float yh[num_points],yh2[num_points]; calculate(M,yh); calculate(2*M,yh2); // doubled step for better accuracy // epsilon is difference between solutions with single // and double steps printf("X\t\tYH\t\tYH2\t\tEPSILON\n"); for (int i=0;i<num_points;i++) printf("%f\t%f\t%f\t%f\n",(x0+((i+1)*(b-a))/num_points), yh[i],yh2[i],(yh2[i]-yh[i])); } -------------------- Demons are a girl's best friends
-------------------- Подарки: (Всего подарков: 171 ) |
|
|
|
|
Текстовая версия | Сейчас: 28.03.2024 - 17:47 |