Архитектура учебник (book.pdf)

Домашнее задание по предмету: «Архитектура компьютеров».

Выполнила студентка 1 курса группы ИУ9-12 Юсупова Лилия.

2014 г.

Постановка задачи

Задача: разработать подпрограмму численного определения длины кривой, заданной в параметрическом виде и, по возможности, улучшить результаты ее работы.

Тип кривой: плоская кривая Безье третьего порядка.

Тип данных: float.

Рисунок 1: График кривой Безье

2

Наивная реализация

Рассмотрим наивную реализацию программы (Листинг 1).

Листинг 1. Наивная реализация.

#include <stdio.h> #include <math.h>

#define reference 2.48112398487 /* эталонное решение */

/* функция вычисления относительной погрешности */ float E (float rez) {

return (rez-reference)/(reference);

}

/* функция вычисления координат */

float F (float a0,float a1,float a2,float a3,float t) { float a;

a=pow((1-t),3)*a0+3*t*pow((1-t),2)*a1+3*t*t*(1-t)*a2+pow(t,3)*a3; return a;

}

/* функция вычисления длины через модуль разности векторов */ float Dl (float X, float x, float Y, float y) {

return (sqrt(pow((X-x),2)+pow((Y-y),2)));

}

/* функция вычисления длины кривой по координатам x,y */

float Dlina (int N,float x0,float y0,float x1,float y1,float x2,float y2,float x3,float y3) { float t,X,Y,x,y,l,L=0;

for (t=0;t<1;t+=(1.0/N)) { x=F(x0,x1,x2,x3,t); y=F(y0,y1,y2,y3,t); l=t+1.0/N; X=F(x0,x1,x2,x3,l); Y=F(y0,y1,y2,y3,l); L+=Dl(X,x,Y,y);

}

return L;

}

void main () {

float x0,y0,x1,y1,x2,y2,x3,y3,L,pogr; int N;

scanf("%f%f",&x0,&y0); /* считывание координат */ scanf("%f%f",&x1,&y1); /* считывание координат */ scanf("%f%f",&x2,&y2); /* считывание координат */ scanf("%f%f",&x3,&y3); /* считывание координат */ printf("Число шагов;Погрешность\n");

/* вывод таблицы для определения зависимости погрешности от кол-ва разбиений */ for (N=1;N<101;N++) { /* цикл для количества разбиений (1..100) */

L=Dlina(N,x0,y0,x1,y1,x2,y2,x3,y3); pogr=E(L); printf("%5d;%20.10f\n", N,pogr);

}

}

3

Рисунок 2: График относительной погрешности для Листинга 1

По графику (Рисунок 2) видно, что при некоторых N (например, N=11) погрешность резко увеличивается.

Гипотеза №1.

Скорее всего так получается потому, что при данных значения N происходит лишняя итерация цикла.

Математическое обоснование:

Попробуем понять, почему же при вычислениях возникают и накапливаются погрешности. Главная причина кроется в устройстве типа float. По стандарту IEEE 754-2008 (включающий в себя стандарт IEEE 754-1985) ширина битового поля, выделенного под экспоненту типа float (равно, как и для типа double) ограничен (Рисунок 3), следовательно, при выполнении арифметических операций происходит округление и, соответственно, мы получаем неточный результат.

На Рисунке 3: S - бит знака (если S=0 – положительное число; S=1 – отрицательное число), E – смещенная экспонента двоичного числа, M – остаток мантиссы двоичного нормализованного числа с плавающей точкой.

Как следствие, получаем, что арифметические операции являются основной причиной возникновения ошибок.

Таким образом, любую арифметическую операцию с типом float можно представить как:

( A B)(1±e) , где e – погрешность операции.

Ожидаемый цикл:

tожид=(((t+dt)+dt)+dt+...+dt )

Полученный цикл:

tитог=((((t+dt )(1±e )+dt)(1±e)+dt)+...+dt )(1±e )

Понятно, что tожид≠tитог , так как e≠0 .

Проверим, так ли это на самом деле.

Настроим отладочную печать (Листинг 2) при N=10 и N=11.

Листинг 2. Настройка отладочной печати.

/* функция вычисления длины кривой по координатам x,y */

float Dlina (int N,float x0,float y0,float x1,float y1,float x2,float y2,float x3,float y3) { float t,X,Y,x,y,l,L=0;

int i=1;

for (t=0;t<1;t+=(1.0/N),i++) { /*добавлена новая переменная i для определения кол-ва шагов*/ x=F(x0,x1,x2,x3,t);

y=F(y0,y1,y2,y3,t); l=t+1.0/N; X=F(x0,x1,x2,x3,l); Y=F(y0,y1,y2,y3,l); L+=Dl(X,x,Y,y);

/*отладочная печать*/

printf("%2d %11.9f %11.9f %11.9f %11.9f %11.9f %11.9f \n", i,x,y,X,Y,t,L);

}

return L;

}

5

Рисунок 4: Отладочная печать для N=10 и N=11

Проанализируем полученный результат.

Как видим из иллюстрации, представленной выше (Рисунок 4), при N=10 программа разбивает кривую на 10 частей, а при N=11 (место «выброса») программа разбивает кривую на 12 частей, а не на 11, как должно быть. Так же, видим, что при N=11 на 11 шаге значение переменной L больше всего похоже на эталонное решение (2.48112398487) .

Возможное решение:

Попробуем решить эту проблему, установив счетчик шагов (Листинг 3).

6

Рисунок 5: Относительная погрешность после установления счетчика

(Листинг 3)

Как можно заметить из графика (Рисунок 5), «выбросы» исчезли. Следовательно, наше решение было верным.

Для того, чтобы график вычисления относительной погрешности был нагляднее, построим его для N=20..100 (Рисунок 6).

7

Рисунок 6: График относительной погрешности

Ограничив количество итераций цикла, мы, тем не менее, не избавились от погрешности вычисления t, а следовательно и не изменили точность

На каждой итерации цикла результат вычисления переменной t суммируется с предыдущим, погрешность накапливается по принципу «снежного кома».

Математическое обоснование:

Каждое последующее слагаемое складывается с предыдущим, который, в свою очередь, уже был получен с погрешностью.

См. Математическое обоснование гипотезы 1, с 4.

Возможное решение:

Заменить суммирование умножением (Листинг 4), на i-том шаге , вместо i суммирований, мы получим всего 1 умножение.

8

Листинг 4. Замена суммирования умножением.

/* функция вычисления длины кривой по координатам x,y */

float Dlina (int N,float x0,float y0,float x1,float y1,float x2,float y2,float x3,float y3) { float t,X,Y,x,y,l,L=0;

int i=0;

for (;t<1 && i<N;i++ /* t+=1.0/N */) { /* убрали суммирование */ t=((float)i)/N; /* изменение переменной t */ x=F(x0,x1,x2,x3,t);

y=F(y0,y1,y2,y3,t);

l=((float)(i+1))/N; /* замена суммирования 2 */

X=F(x0,x1,x2,x3,l);

Y=F(y0,y1,y2,y3,l); L+=Dl(X,x,Y,y);

}

return L;

}

Теперь мы не накапливаем сумму при каждой итерации цикла, а вычисляем значение переменной t, пользуясь всего одной арифметической операцией.

Рисунок 7: Относительная погрешность до и после изменения программы.

Оценим результат внесенных изменений. Судя по графикам (Рисунок 7), модуль относительной погрешности уменьшился. Следовательно, решение заменить суммирование умножением было правильным.

9

Рисунок 8: Результаты работы программы при N=100..5000

Попробуем изменить начальные условия так, чтобы этот эффект («скачки») проявлялся при малых N. Подберем начальные условия такие, чтобы отношение длин минимального и максимального вектора при заданном разбиении было

равно некоторому коэффициенту k. Для этого закрепим три вектораP0 , P1 , P3 и найдем их длину:

P0=√0+0=0 ,

P1=√1+1=√2 ,

P3=√1+0=1 .

Вектор P2 оставим в виде коэффициента. Запишем полученную формулу для вычисления длины вектора r :

10