Шпаргалка по языку СИ. Конспект Лекций «программирование На Языке Высокого Уровня Си» П. Конспект лекций для студентов высших учебных заведений, обучающихся по специальности 230102. 65 Автоматизированные системы обработки информации и управления
 Скачать 1.25 Mb.
|
|
{ int n = 4; printf("Введите количество дисков : "); scanf("%d",&n); printf("\nХаной %d дисками :\n", n); hanoy(1, 3, n, 0); getch(); // system ("pause"); } Результат выполнения программы Введите количество дисков : 4 Ханой с 4 дисками : перемещение(1,2) перемещение(1,3) перемещение(2,3) перемещение(1,2) перемещение(3,1) перемещение(3,2) перемещение(1,2) перемещение(1,3) перемещение(2,3) перемещение(2,1) перемещение(3,1) перемещение(2,3) перемещение(1,2) перемещение(1,3) перемещение(2,3) Пример 13. Функция решение кубического уравнения с действительными коэффициентами методом Виета-Кардано. Корни могут быть комплексными. Кубическое уравнение записывается в виде: x 3 +a*x 2 +b*x+c=0. Для нахождения его корней, в случае действительных коэффициентов, вначале вычисляются: Q = (a2-3b)/9 R = (2a3-9ab+27c)/54. Далее, если R2 < Q3, то уравнение имеет три действительных корня, вычисляющихся по формулам (Виета): е = acos(R/sqrt(Q3))/3, x 1 = -2*sqrt(Q)cos(t)-a/3, x 2 = -2*sqrt(Q)cos(t+(2*pi/3))-a/3, x 3 = -2*sqrt(Q)cos(t-(2*pi/3))-a/3. В том случае, когда R2 >= Q3, то действительных корней один (общий случай) или два (вырожденные случаи). Кроме действительного корня, имеется два комплексно-сопряженных. Для их нахождения вычисляются (формула Кардано): 202 A = -sign(R)[|R|+sqrt(R2-Q3)]1/3, B = Q/A при A != 0 или B = 0 при A = 0. Действительный корень будет: x 1 = (A+B)-a/3. Комплексно-сопряженные корни: x 2,3 = -(A+B)/2-a/3 + i*sqrt(3)*(A-B)/2 В том случае, когда A = B, то комплексно-сопряженные корни вырождаются в действительный: x2=-A-a/3. Формулы Кардано и Виета требуют применения специальных функций, и в том случае, когда требуется провести большую серию вычислений корней кубического уравнения с не слишком сильно меняющимися коэффициентами, более быстрым алгоритмом является использование метода Ньютона или других итерационных методов (с нахождением начального приближения по формулам Кардано-Виета). #include #include #include #include #include #define M_PI (3.141592653589793) #define M_2PI (2.*M_PI) /* int Cubic(double *x,double a,double b,double c); Параметры: a, b, c - коэффициенты x - массив решения (size 3). На выходе: 3 действительных корня -> затем x ими заполняется; 1 действительный + 2 комплексных -> x[0] - действительный, x[1] действительная часть комплексных корней, x[2] - неотрицательная мнимая часть. Результаты: 1 - 1 действительный + 2 комплексных; 2 - 1 действительный корень + мнимая часть комплексных корней, если 0 (т.е. 2 действительных корня). 3 - 3 действительных корня; */ int Cubic (double *x, double a, double b, double c) { double q, r, r2, q3; q = (a * a - 3. * b) / 9.; r = (a * (2. * a * a - 9. * b) + 27.*c) / 54.; r2 = r * r; q3 = q * q * q; if (r2 < q3) { double t = acos(r / sqrt (q3)); 203 a /= 3.; q = -2. * sqrt (q); x[0] = q * cos (t / 3.) - a; x[1] = q * cos ((t + M_2PI) / 3.) - a; x[2] = q * cos ((t - M_2PI) / 3.) - a; return (3); } else { double aa, bb; if (r <= 0.) r =- r; aa = -pow (r + sqrt (r2 - q3), 1. / 3.); if(aa != 0.) bb = q / aa; else bb = 0.; a /= 3.; q = aa + bb; r = aa - bb; x[0] = q - a; x[1] = (-0.5) * q - a; x[2] = (sqrt (3.) * 0.5) * fabs (r); if (x[2] == 0.) return (2); return (1); } } Пример 14. Алгоритм Краута (нижне-верхняя (LU) декомпозиция матрицы). Алгоритм Краута − это фактически метод решения систем общего вида, конкурирующий по быстродействию с общеизвестным методом Гаусса- Жордана, но позволяющий более эффективно использовать решение. Если мы можем разложить матрицу линейной системы A в произведение A=L*U, где L − нижняя, а U − верхняя треугольные матрицы, то решение системы уравнений с произвольной правой частью производится весьма просто, применением двух обратных подстановок. Более того, в отличие от известного метода Гаусса-Жордана, разложенная матрица позволяет быстро решать серии линейных уравнений с различными правыми частями при одной и той же матрице. Метод Краута позволяет провести LU-декомпозицию матрицы примерно за то же число операций, что и "прямая" часть метода Гаусса-Жордана. Итоговые коэффициенты двух треугольных матриц упаковываются в матрицу того же размера, что и A, и на том же месте в памяти; при этом верхняя матрица U размещается в наддиагональной части и на диагонали; нижняя L в 204 поддиагональной части, а диагональные элементы L считаются все равными 1 (без ограничения общности) и не выводятся. Алгоритм Краута представляет из себя следующее: 1. Положить L ii =1 i=0,...,N-1 (здесь ничего не производится, а только имеется в виду для дальнейших шагов; 2. Для каждого j=0,1,...,N-1 проделать: 1. Для всех i=0,...,j вычислить Uij: U ij =A ij - SUM(0<=k ik *U ik ) (при i=0 сумму полагать нулем); 2. Для всех i=j+1,...,N-1 вычислить: L ij = A ij - SUM(0<=k ik *U ik ))/U ii Обе процедуры выполняются до перехода к новому j. Устойчивая работа алгоритма Краута возможно только при условии выбора ведущего элемента для каждого обращения к j из п.2 алгоритма. Выбор производится перед выполнением п.2 для очередного j перестановки необработанных строк матрицы так, чтобы строка, подставляемая на место j, содержала наибольший по модулю элемент в колонке j. Порядок перестановки необходимо запомнить, для чего достаточно использовать дополнительный вектор целых величин. Еще один вектор используется внутри алгоритма для масштабирования. Алгоритм Краута имеет несколько приложений, одно из которых − решение системы линейных уравнений (обратная подстановка с учетом порядка перестановки строк), другие − вычисление обратной матрицы и вычисление детерминанта. Подробнее вызовы функций, связанных с алгоритмом Краута и реализация деталей алгоритма находятся в комментариях к программе и в ней самой. #include #include #include #include #include #define TINY 1.e-30 /* LU-декомпозиция в соответствии с алгоритмом Краута. Описание: 205 int LU_decompos(double **a,int n,int *indx,int *d,double *vv); Параметры: a - исходная матрица (n x n) на входе; n - размер матрицы; indx - массив целых чисел (размера n) для запоминания перемещений; d - на выходе, содержит +1 или -1 для четного или нечетного числа перемещений. vv - временный массив (size n). Результат: 0 - некорректная входная матрица (непригодна для декомпозиции), 1 - положительный результат. Обратная замена, использующая LU-декомпозицию матрицы. Описание: void LU_backsub(double **a,int n,int *indx,double *b); Параметры: a - матрица, разложенная по Крауту; n - размер матрицы; indx - порядок перемещения, полученный алгоритмом декомпозиции; b - вектор (размера n). Примечание: a и indx не изменяются в этой функции и могут быть использованы повторно. Инвертирование матрицы с использованием LU-разложенной матрицы. Описание: void LU_invert(double **a,int n,int *indx,double **inv,double *col); Параметры: a - матрица, разложенная по Крауту; n - размер матрицы; indx - порядок перестановки, полученный алгоритмом декомпозиции; inv - матрица назначения; col - временный массив (размера n). Получение матрицы, полученной с использованием LU-разложенной матрицы Описание: double LU_determ(double **a,int n,int *indx,int *d); Параметры: a - матрица, разложенная по Крауту; n - размер матрицы; indx - порядок перемещения, полученный алгоритмом декомпозиции; d - знак равенства (+1 или -1) полученный при декомпозиции. Результат: определяющее значение. */ /* декомпозиция */ int LU_decompos (double **a, int n, int *indx, int *d, double *vv) { register int i, imax, j, k; double big, sum, temp; for (i = 0; i < n; i++) { big = 0.; for (j = 0; j < n; j++) if ((temp = fabs (a[i][j])) > big) big = temp; 206 if (big == 0.) return (0); vv[i] = big; } /* главный цикл алгоритма Краута */ for (j = 0; j < n; j++) { /* это часть a) алгоритма, исключая i==j */ for (i = 0; i < j; i++) { sum = a[i][j]; for (k = 0; k < i; k++) sum -= a[i][k] * a[k][j]; a[i][j] = sum; } /* инициализация для поиска наибольшего элемента */ big = 0.; imax = j; for (i = j; i < n; i++) { sum = a[i][j]; for (k = 0; k < j; k++) sum -= a[i][k] * a[k][j]; a[i][j] = sum; if ((temp = vv[i] * fabs (sum)) >= big) { big = temp; imax = i; } } /* обмен строк, если нужен, смена равенства и размера множителя */ if (imax != j) { for (k = 0; k < n; k++) { temp = a[imax][k]; a[imax][k] = a[j][k]; a[j][k] = temp; } *d =- (*d); vv[imax] = vv[j]; } /* сохраняем индекс */ indx[j] = imax; if (a[j][j] == 0.) a[j][j] = TINY; if (j < n - 1) { temp = 1. / a[j][j]; for (i = j+1; i < n; i++) a[i][j] *= temp; } } return (1); } /* обратная замена */ void LU_backsub (double **a, int n, int *indx, double *b) { register int i, j, ip, ii = -1; double sum; /* первый шаг обратной замены */ for (i = 0; i < n; i++) { ip = indx[i]; sum = b[ip]; b[ip] = b[i]; if (ii >= 0) for (j = ii; j < i; j++) sum -= a[i][j] * b[j]; 207 else if (sum) ii = i; /* получен ненулевой элемент */ b[i] = sum; } /* второй шаг */ for (i = n - 1; i >= 0; i--) { sum = b[i]; for (j = i + 1; j < n; j++) sum -= a[i][j] * b[j]; b[i] = sum / a[i][i]; } } /* инвертирование матрицы */ void LU_invert (double **a, int n, int *indx, double **inv, double *col) { register int i, j; for (j = 0; j < n; j++) { for (i = 0; i < n; i++) col[i] = 0.; col[j] = 1.; LU_backsub (a, n, indx, col); for (i = 0; i < n; i++) inv[i][j] = col[i]; } } /* определяющее вычисление*/ double LU_determ (double **a, int n, int *indx, int *d) { register int j; double res = (double)(*d); for (j = 0; j < n; j++) res *= a[j][j]; return (res); } Пример 15. Алгоритм поиска интервала нахождения корня нелинейной функции. Выяснение интервала, на котором корень содержится является важной проблемой поиска корня нелинейной функции действительной переменной. Здесь приведен алгоритм поиска такого интервала и ограничения на его применение. Примем, что корень функции f(x) окружен на интервале [a,b], если f(a) и f(b) имеют противоположные знаки. Чтобы окруженный согласно этому определению корень действительно существовал на этом интервале, достаточно непрерывности f(x), а для его единственности − еще и монотонности. При невыполнении этих свойств возможно отсутствие корня на [a,b] или неопределенность его позиции. При использовании компьютера, всегда имеем дело с дискретным набором возможных представлений чисел (хотя и достаточно плотным). Кроме того, монотонность вычисленной функции может быть слегка нарушена в пределах 208 точности ее вычисления. Это в ряде случаев усложняет вычисление окруженных корней функции, если к их точности предъявляются завышенные требования. Окружение корня функции при гарантии ее определения на неограниченном интервале, производится по следующему итерационному алгоритму. 1. Для начального приближения x 0 , найти f 0 =f(x 0 ), задать начальный интервал поиска D и его инкремент d>1. 2. Вычислить a=x 0 -D, b=x 0 +D; f a =f(a), f b =f(b). 3. Увеличить интервал поиска: D=D*d. Если интервал превысил некоторый заданный предел, то выйти с индикацией ошибки. 3.a. Если знаки f a и f 0 отличаются, то считать корень окруженным на [a,x 0 ] − > выход. 3.b. Если знаки f b и f 0 отличаются, то считать корень окруженным на [x 0 ,b] -> выход. 3.c. Если f 0 > 0 (случай меньше нуля делается аналогично) алгоритм продолжается: 4. Проверяется, какое из f a или f b наименьшее. Если оба одинаковы, то переходим к 4a (двусторонний поиск), если f b − производим поиск вправо 4b, иначе − поиск влево 4c. 4.a. Находим a=a-D, b=b+D, f a =f(a), f b =f(b), идем к пункту 3. 4.b. Присваиваем последовательно a=x 0 , x 0 =b, f a =f 0 , f 0 =f b ; находим b=b+D, f b =f(b), идем к пункту 3. 4.c. Аналогично 4b, только направление поиска − влево. Так как интервал поиска постоянно расширяется, то в конце концов используя указанный алгоритм корень будет окружен. Возможны модификации алгоритма в двух направлениях: 1. Увеличивать интервал не в геометрической прогрессии, а в арифметической либо по заданному сценарию; 209 2. Если область определения функции заведомо ограничена, то расширение интервала поиска также следует ограничивать имеющимися пределами, либо доопределять функцию там, где ее оригинал не определен. #include #include #include #include #include /* Предполагается, что функция существует на бесконечном и непрерывном промежутке. int BracketRoot (double x0, double *a, double *b, double d0, double di, double dmax, double (*fun)(double)); Параметры: x0 - начальное приближение; a - левая граница; b - правая граница; d0 - начальный интервал поиска; di - интервал расширения (в геометрической прогрессии); dmax - максимальный интервал; fun - указатель на функцию. Возвращает: 1 - если корень окружен; 0 - сбой */ int BracketRoot (double x0, double *a, double *b, double d0, double di, double dmax, double (*fun)(double)) { double fa, fb, f0; f0 = (*fun)(x0); *a = x0 - d0; *b = x0 + d0; fa = (*fun)(*a); fb = (*fun)(*b); while ((d0 *= di) < dmax) { if (f0 >= 0.) { if (fa < 0.) { *b = x0; return (1); } if (fb < 0.) { *a = x0; return (1); } if (fa > fb) { *a = x0; x0 = (*b); *b += d0; fa = f0; f0 = fb; 210 fb = (*fun)(*b); } else if (fa < fb) { *b = x0; x0 = (*a); *a -= d0; fb = f0; f0 = fa; fa = (*fun)(*a); } else { *a -= d0; *b += d0; fa = (*fun)(*a); fb = (*fun)(*b); } } else if (f0 < 0.) { if (fa >= 0.) { *b = x0; return (1); } else if (fb >= 0.) { *a = x0; return (1); } if (fa < fb) { *a = x0; x0 = (*b); *b += d0; fa = f0; f0 = fb; fb = (*fun)(*b); } else if (fa > fb) { *b = x0; x0 = (*a); *a -= d0; fb = f0; f0 = fa; fa = (*fun)(*a); } else { *a -= d0; *b += d0; fa = (*fun)(*a); fb = (*fun)(*b); } } } return (0); } 211 Пример 16. Алгоритм нахождения Полинома Лагранжа (интерполяция Лагранжа). Интерполяционный многочлеен Лагранжа − многочлен минимальной степени, принимающий данные значения в данном наборе точек. Для n + 1 пар чисел (x_0, y_0), (x_1, y_1) \ dots (x_n, y_n), где все x i различны, существует единственный многочлен L(x) степени не более n, для которого L(x i ) = y i В простейшем случае (n = 1) − это линейный многочлен, график которого − прямая, проходящая через две заданные точки. |