16237

Визуализация численных методов. Решение обыкновенных дифференциальных уравнений

Курсовая

Информатика, кибернетика и программирование

Визуализация численных методов. Решение обыкновенных дифференциальных уравнений. Содержание Введение. 1. Постановка задачи и математическая модель. 2. Описание численных методов применительно к конкретной задаче 3. Блоксхемы программ и основных подпрограм

Русский

2013-06-20

350 KB

5 чел.

Визуализация численных методов.

Решение обыкновенных дифференциальных уравнений.

Содержание

Введение. 

1. Постановка задачи и математическая модель. 

2. Описание численных методов (применительно к конкретной задаче) 

3. Блок-схемы программ и основных подпрограмм. 

4. Листинг программы на языке VisualBasic. 

5. Формы проекта

6. Решение задачи в Mahtcad. 

Заключение. 


Введение.

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

Дифференциальным называется уравнение, содержащее одну или несколько производных. Лишь очень немногие из них удается решить без помощи вычислительной техники. Поэтому численные методы решения дифференциальных уравнений играют огромную роль в практике инженерных расчетов.

Для решения таких уравнений удобно использовать языки программирования, так как они позволяют быстро и точно найти решения уравнения и построить график интегральной кривой.

Целью данной работы является решение задачи Коши для дифференциального уравнения 1-го порядка.

Чтобы достичь поставленной цели были выдвинуты следующие задачи:

1. Изучить численные методы решения дифференциальных уравнений;

2. Самостоятельно вычислить первую точку интегральной кривой заданными методами;

3. Написать программу для построения интегральной кривой на языке программирования VusuaiBasic;

4. Проверить решение в среде MathCad.

1. Постановка задачи и математическая модель.

Решить методами Эйлера и Эйлера модифицированного задачу Коши для дифференциального уравнения 1-го порядка на отрезке [X0; Xk] с шагом h и начальным условием: Y(X0) = Y0.

Ответ должен быть получен в виде таблицы результатов:

X

Y(1)

Y(2)

YT

X0

Y0(1)

Y0(2)

Y(X0)

X1

Y1(1)

Y1(2)

Y(X1)

Xk

Yk(1)

Yk(2)

Y(Xk)

Где Y(1), Y(2) – решения, полученные различными численными методами, YT – точное решение дифференциального уравнения.

Возможно представление результатов решения не в виде таблицы, а в виде списков.

Данные таблицы визуализировать на форме в виде графиков.

Перед вычислением последнего столбца таблицы результатов необходимо из начальных условий вычесть значение коэффициента c, используемого в общем решении.

  

Дифференциальное уравнение

X0

Xk

h

Y0

Общее решение

(y^2-2*x*y)dx+x^2=0

1

2

0,1

0.2

y=x^2/(c+x)


2. Описание численных методов (применительно к конкретной задаче)

Метод Эйлера модифицированный

Этот метод часто используют для уменьшения погрешности вычислений.

Пусть дано дифференциальное уравнение первого порядка

y’= f(x,y)

с начальным условием

y(x0)=y0

выберем шаг h и введем обозначения:

xi = x0 +i h и yi = y(xi ), где i=0,1,2… 

xi  узлы сетки,

yi  значения интегральной функции в узлах

При использовании модифицированного метода Эйлера шаг h делится на два отрезка.

Метод Эйлера модифицированный


Расчет первой точки методом Эйлера модифицированного.

  1.  строим оси координат;
  2.  отмечаем начальную точку интегральной кривой А(x0;y0);
  3.  ищем угол наклона касательной к графику в точке А(1;0.2);

y’(x0) = f(x0;y0)

α0 = arctg(f(x0;y0))

α0 = arctg((0.2/1)*( 2-0.2/1)) = arctg(0.36)

  1.  строим касательную L0 к графику функции в точке А под углом α0 ;
  2.  находим x1 по формуле x1=x0+i*h

x1=1+1*0.1=1.1

  1.  отмечаем середину отрезка x1/2=x0+h/2 и проводим прямую из этой точки до прямой L0, отмечаем точку B(xb;yb);
  2.  ищем координаты точки В:

xb = x0+h/2

xb = 1+0.05=1.05

yb =  y0+h/2*f(x0;y0)

yb = 2+0.05*((0.2/1)*(2-0.2/1))

yb = 2.018

значит точка В имеет координаты (1.05;2.018);

  1.  ищем угол наклона касательной к графику в точке В:

αb = arctg(f(xb;yb))

αb= arctg((2.018/1.1)*(2-2.018/1.1)=arctg(0.304);

  1.  строим касательную L1 в точке В под углом αb;
  2.  проводим прямую x1=x0+h (=1+0.1=1.1)до пересечения с прямой L1, отмечаем точку С(x1;y1);
  3.  ищем y точки С:

y1=yb+h/2*f (xb;yb)

y1= 2.018+0.05*0.304

y1=2.17

значит точка С имеет координаты (1.1; 2.17)

Метод Рунге-Кутта 4-го порядка

Пусть дано дифференциальное уравнение первого порядка

y’= f(x,y)

c начальным условием

y(x0)=y0

выберем шаг h и введем обозначения:

xi = x0 +i h и yi = y(xi ), где i=0,1,2… 

Согласно методу Рунге-Кутта четвертого порядка, последовательность значения yi  искомой функции y определяется по формуле:

yi+1 = yi +∆ yi

где

∆ yi= 1/6* ( k1+2k2+2k3+k4 ), i=0,1,2…

а числа k1 ω , k2 ω, k3 ω , k4 ω на каждом шаге вычисляются по формулам:

k1= h * f ( x1, y1 )

k2= h * f ( x1+h/2, y1+ k1/2 )

k3= h * f ( x1+h/2, y1+ k2/2 )

k4= h * f ( x1+h, y1+k3 )

Метод Рунге-Кутта легко программируется и обладает значительной точностью и устойчивостью для широкого круга задач.

Метод Рунге-Кутта 4 порядка
Расчет первой точки методом Рунге-Кутта 4 порядка.

k1 =  h*f(x0;y0)

k1 = 0.1*(0.2*(2-0.2)) = 0.036

k2 = h*f(x0+h/2; y0+k1/2)

k2 = 0.1*((0.218/1.05)*(2-0.218/1.05)) = 0.037

k3 = h*f(x0+h/2;y0+k2/2)

k3 = 0.1*((0.2185/1.05)(2-0.2185)) = 0.0373

k4 = h*f(x0+h/2;y0+k3)

k4 = 0.1*((0.2373/1.1)(2-0.2373/1.1))=0.038

∆y = k = (k1+2*k2+2*k3+k4)/6

y = k = (0.036+2*0.037+2*0.0373+0.038)/6 = 0.036

y1=∆y +y0

y1 = 0.036 + 0.2 = 0.236

3. Блок-схемы программ и основных подпрограмм.

Для решения задачи и составления программы составим основные блок-схемы, которые приведены ниже.

1.Подпрограмма метода Эйлера модифицированного

 

 


2. Подпрограмма метода Рунге-Кутта 4 порядка

3. Подпрограмма общего решения функции

4. Алгоритм функции


5. Алгоритм программы



4. Листинг программы на языке VisualBasic.

Dim x(), em(), rk(), o() As Single

Private i, n As Integer

Private x0, xk, y0, h, miny, maxy, minx, maxx As Single

Function f(a, b) As Single

f = -((b * b) - 2 * a * b) / (a * a)

End Function

Private Sub EilerM()

ReDim x(n)

ReDim em(n)

x(0) = x0

em(0) = y0

For i = 0 To n - 1

x(i) = Round(x0 + i * h, 3)

em(i + 1) = Round(em(i) + h * f(x(i) + h / 2, em(i) + h / 2 * f(x(i), em(i))), 3)

Next i

End Sub

Private Sub RungeKutt()

ReDim x(n)

ReDim rk(n)

x(0) = x0

rk(0) = y0

For i = 0 To n - 1

x(i) = Round(x0 + i * h, 3)

 k1 = h * f(x(i), rk(i))

 k2 = h * f(x(i) + (h / 2), rk(i) + (k1 / 2))

 k3 = h * f(x(i) + (h / 2), rk(i) + (k2 / 2))

 k4 = h * f(x(i) + h, rk(i) + k3)

 k = (k1 + 2 * k2 + 2 * k3 + k4) / 6

rk(i + 1) = Round(rk(i) + k, 3)

Next i

End Sub

Private Sub Obchee()

ReDim x(n)

ReDim o(n)

c = Round(((x0 * x0 / y0) - x0), 3)

For i = 0 To n

x(i) = Round(x0 + (i * h), 3)

o(i) = Round((x(i) * x(i)) / (c + x(i)), 3)

Next i

End Sub

Private Sub Command1_Click()

x0 = Val(Text1.Text)

xk = Val(Text2.Text)

y0 = Val(Text3.Text)

h = Val(Text4.Text)

n = Round((xk - x0) / h)

c = Round(((x0 * x0 / y0) - x0), 3)

MSFlexGrid1.Cols = 4

MSFlexGrid1.Rows = n + 2

MSFlexGrid1.TextMatrix(0, 0) = "x"

MSFlexGrid1.TextMatrix(0, 1) = "Obchee"

MSFlexGrid1.TextMatrix(0, 2) = "EilerM"

MSFlexGrid1.TextMatrix(0, 3) = "RungeKutt"

EilerM

RungeKutt

Obchee

For i = 0 To n

MSFlexGrid1.TextMatrix(i + 1, 0) = Str(x(i))

MSFlexGrid1.TextMatrix(i + 1, 1) = Str(o(i))

MSFlexGrid1.TextMatrix(i + 1, 2) = Str(em(i))

MSFlexGrid1.TextMatrix(i + 1, 3) = Str(rk(i))

Next i

miny = o(0)

maxy = o(n)

For i = 1 To n

If em(i) > maxy Then maxy = em(i)

If em(i) < miny Then miny = em(i)

If rk(i) > maxy Then maxy = rk(i)

If rk(i) < miny Then miny = rk(i)

If o(i) > maxy Then maxy = o(i)

If o(i) < miny Then miny = o(i)

Next i

Label9.Caption = Str(maxy)

Label10.Caption = Str(maxx)

Label11.Caption = Str(miny)

Label12.Caption = Str(minx)

Picture1.Cls

kx = (3600 - 240) / (xk - x0)

ky = (5640 - 360) / (maxy - miny)

For i = 0 To n - 1

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (em(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (em(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(0, 0, 9999)

Next i

For i = 0 To n - 1

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (rk(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (rk(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(0, 9999, 0)

Next i

For i = 0 To n - 1

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (o(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (o(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(9999, 0, 0)

Next i

End Sub

Private Sub Command2_Click()

End

End Sub

5. Формы проекта

Форма разработки проекта.

 


Форма с результатами работы проекта.


6. Решение задачи в Mathcad.


 


Заключение.

В курсовой работе были описаны методы Эйлера модифицированного и Рунге-Кутта четвертого порядка. С помощью которых была решена задача Коши первого порядка. В ходе работы были построены блок-схемы программ, разработан код программы построения графика в VisualBasic.

С помощью данной работы были закреплены навыки работы в различных приложениях WINDOWS, а именно Microsoft Word, Mathcad и  VisualBasic.

Курсовая работа на тему «Визуализация численных методов» показалась мне интересной и была выполнена по указаниям преподавателя.


y

x

y=y(x)

О

А

С

В

h/2

h

xi

xi+1

ε1

ε

α

α1

EilerM(x0,xk,y0,n,em)

h=(xk-x0)/n

i = 1, …, n-1

x = x0 + i * h

em(i + 1) = em(i) + h * f(x + h / 2, em(i) + h / 2 * f(x, em(i))

End

RungeKutt (x0,xk,y0,n,rk)

h=(xk-x0)/n

i = 1, …, n-1

nd

x(i) = x0 + i * h

k1 = h * f(x, rk )

k2 = h * f(x + (h / 2), rk  + (k1 / 2))

k3 = h * f(x + (h / 2), rk + (k2 / 2))

k4 = h * f(x + h, rk  + k3)

k = (k1 + 2 * k2 + 2 * k3 + k4) / 6

rk (i + 1) = rk (i) + k

Obchee

x(n), o(n)

i = 1, …, n

x(i) = (x0 + (i * h)

o(i) = ((x(i) * x(i)) / (c + x(i))

End

f (a,b)

f=b/a*(2-b/a)

End

Start

x0,xk,y0,h

n=(xk-x0)/h

C=(x^2/y)-x

i=0,..,n-1

x=x0+i*h

em(i)>maxy

Maxy=em(i)

em(i)<miny

Miny=em(i)

да

нет

да

нет

EilerM

RungeKutt

Obchee

rk(i)>maxy

Maxy = rk(i)

rk(i)<miny

miny=rk(i)

o(i)>maxy

maxy=o(i)

o(i)<miny

miny=o(i)

да

нет

да

нет

да

нет

да

нет

i=0,..,n-1

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (em(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (em(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(0, 0, 9999)

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (rk(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (rk(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(0, 9999, 0)

z1 = Round(240 + (x(i) - x0) * kx)

z2 = Round(5640 - (o(i) - miny) * ky)

z3 = Round(240 + (x(i + 1) - x0) * kx)

z4 = Round(5640 - (o(i + 1) - miny) * ky)

Picture1.Line (z1, z2)-(z3, z4), RGB(9999, 0, 0)

End

i=0,..,n-1

MSFlexGride

Label7

Label 6

Label 1

Text3

Text2

Text1

Picture1

Command2

Command1

Text4

Label 5

Label 4

Label 3

Label 2


 

А также другие работы, которые могут Вас заинтересовать

10457. Алгоритмы сжатия на основе вейвлет-преобразования. Алгоритм SPIHT 63 KB
  Алгоритмы сжатия на основе вейвлетпреобразования. Алгоритм SPIHT. Изображение полученное при помощи вейвлетпреобразования можно сжимать различными способами. Большинство из них можно отнести к одной из двух категорий. К первой категории относятся способы сводящиеся
10458. Алгоритмы обработки изображений в астроориентации 666.5 KB
  Алгоритмы обработки изображений в астроориентации. Введение К приборам астроориентации космических аппаратов относятся солнечные датчики датчики положения Земли и звездные датчики. Термин датчик не должен вводить в заблуждение солнечные датчики и датчики полож
10459. Улучшение изображений. Изменение контраста, Видоизменение гистограммы. Подавление шумов. Медианная фильтрация 247.5 KB
  Улучшение изображений. Изменение контраста Видоизменение гистограммы. Подавление шумов. Медианная фильтрация. Процедура улучшения изображений сводится к выполнению комплекса операций с целью либо улучшения визуального восприятия либо преобразования в форму боле...
10461. Операционные системы - основные понятия 79.05 KB
  Операционные системы основные понятия. Операционная система сокр. ОС англ. operating system комплекс управляющих и обрабатывающих программ которые с одной стороны выступают как интерфейс между устройствами вычислительной системы и прикладными программами а с друг
10462. Ядро операционной системы 35.56 KB
  Тема: Ядро операционной системы. Ядро центральная часть операционной системы ОС обеспечивающая приложениям координированный доступ к ресурсам компьютера таким как процессорное время память и внешнее аппаратное обеспечение. Также обычно ядро предоставляет сервис...
10463. Операционные системы реального времени. Архитектуры ОСРВ 56.33 KB
  Тема: Операционные системы реального времени. Операционная система реального времени ОСРВ англ. RealTime Operating System тип операционной системы. Есть много определений термина по сути похожих друг на друга. Самые распространённые из них: Операционная система в ...
10464. Стандарты ОСРВ 37.03 KB
  Тема: Стандарты ОСРВ. Большие различия в спецификациях ОСРВ и огромное количество существующих микроконтроллеров выдвигают на передний план проблему стандартизации в области систем реального времени. Наиболее ранним и распространенным стандартом ОСРВ является...
10465. Настраиваемость операционных систем 69.04 KB
  Тема: Настраиваемость операционных систем. В последнее время одной из главных тем исследовательских работ в области операционных систем стало исследование настраиваемости customizability или адаптируемости операционной системы. Настраиваемой или адаптируемой операци