Численное моделирование и анализ переходных процессов в электрической цепи

Нижегородский государственный технический университет им Р.Е. Алексеева

Курсовая работа по информатике

" Численное моделирование и анализ переходных процессов в электрической цепи"

Выполнил: Воронкин И.С.

группы 11-Э-3

Проверила: Кулагина Л.В.

г. Нижний Новгород — 2012 г.

• 1. Постановка задачи
• 2. Вывод системы дифференциальных уравнений
• 3. Теоретическая часть
• 4. Численное интегрирование
• 4.1 Метод левых прямоугольников
• 4.2 Метод средних прямоугольников
• 4.3 Формула средних прямоугольников
• 4.4 Метод правых прямоугольников
• 4.5 Формула правых прямоугольников
• 4.6 Метод Симпсона
• 4.7 Метод трапеций
• 5. Постановка задачи Коши
• 6. Разностные схемы Эйлера
• 7. Метод Рунге-Кутта второго порядка (Метод Эйлера-Коши)
• 8. Метод Рунге-Кутта четвертого порядка
• 9. Практическая часть
• Выводы
• Список литературы

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

Дана схема электрической цепи, содержащая источник переменного тока, катушку индуктивности, конденсатор, набор резисторов и ключ (рис.1).

Рис. 1.

Параметры элементов цепи:

— гармонический источник тока; = 15 В — амплитуда колебаний; - циклическая частота; f, Гц — линейная частота; - фаза; t — текущее время; = 30 Ом, = 25 Ом, = 50 Ом, = 1,88 Ом, = 15 Ом, = 50 Ом — резисторы; L = 5,57 мГн — катушка индуктивности; C = 20 мкФ — конденсатор. Параметры f, для данного варианта принимают следующие значения: f = 40 Гц;. =4п/5

В начальный момент времени ключ находится в положении 1. При этом цепь разомкнута, напряжение на конденсаторе и ток в катушке равны нулю (U = 0, I = 0). Происходит первое переключение ключа (ключ мгновенно переводится в положение 2). При этом происходит заряд конденсатора, меняются значения U и I.

В момент времени ключ мгновенно переключается в положение 1. Конденсатор разряжается, вновь меняются значения U и I. Анализ схемы заканчивается в момент .

2. Вывод системы дифференциальных уравнений

В соответствии с рисунком запишем выражения для I и II законов Кирхгофа для положения ключа 1:

(1)

Систему (1) можно преобразовать, исключив токи I₁ и I₂. Тогда для величин I и U получим систему дифференциальных уравнений первого порядка:

(2)

Аналогично может быть получена система дифференциальных уравнений для величин I и U при положении ключа 2. В этом случае имеем:

(3)

В интервале решается система (3) с начальными условиями:

; В интервале решается система (2). В качестве начальных условий для системы (2), следует использовать соответствующие значения, полученные в результате решения системы (3).

3. Теоретическая часть

1. Аппроксимация — это задача, в результате решения которой находят некоторую аппроксимирующую функцию f (х), такую, чтобы отклонения ее от заданной табличной функции было наименьшим. Чаще всего функцию f (х) представляют в виде полинома по степеням х. Общий вид полинома n-ой степени: f (x) =a₀+a₁x+a₂x²+…+a_nxⁿ.

2. Метод наименьших квадратов. Пусть общее количество точек равно m. Неизвестные коэффициенты а₀, а₁,…a_n, n находим из условия минимизации суммы квадратов отклонений искомой функции от исходных точек. Опуская промежуточные преобразования получим систему уравнений: ZЕA=B, где Z — квадратная матрица размерностью (n+1) x (n+1), составленная из известных координат точек, А — вектор неизвестных коэффициентов; В — вектор-столбец свободных членов (i=1,m).

;; (1)

3. Интерполяция — является частным случаем аппроксимации. Это задача о нахождении такой аналитической функции f (х), которая принимает в точках (узлах) x_i заданные значения y_i

4. Метод неопределенных коэффициентов

Пусть табличная функция содержит m точек. В этом случае можно построить различные виды кусочной интерполяции (кусочно-линейная, кусочно-параболическая и т. д.). В случае непрерывной интерполяции, когда используются все точки одновременно, функцию f (х) будем искать в виде полинома степени n: f (x) =a₀+a₁x+a₂x²+…a_nxⁿ. Степень полинома всегда на единицу меньше числа точек. Следовательно, справедливо соотношение: n=m-1.

Для нахождения неизвестных коэффициентов необходимо построить систему линейных уравнений m-го порядка из условия прохождения полинома через все m точек:

Данную систему можно решить методом Гаусса, Зейделя, Простой интерации, а также использованием специальных утилит.

;; (2)

4. Численное интегрирование

4.1 Метод левых прямоугольников

Разделим отрезок [a; b] на n равных частей, т. е. на n элементарных отрезков. Длина каждого элементарного отрезка. Точки деления будут: x₀=a; x₁=a+h; x₂=a+2Ч h,., x_n-1=a+ (n-1) Ч h; x_n=b. Эти числа будем называть узлами. Вычислим значения функции f (x) в узлах, обозначим их y₀, y₁, y₂,., y_n. Cталобыть, y₀=f (a), y₁=f (x₁), y₂=f (x₂),., y_n=f (b). Числа y₀, y₁, y₂,., y_n являются ординатами точек графика функции, соответствующих абсциссам x₀, x₁, x₂,., x_n. Площадь криволинейной трапеции приближенно заменяется площадью многоугольника, составленного из n прямоугольников. Таким образом, вычисление определенного интеграла сводится к нахождению суммы n элементарных прямоугольников.

Формула левых прямоугольников:

Рис. 1.

4.2 Метод средних прямоугольников

Разделим отрезок [a; b] на n равных частей, т. е. на n элементарных отрезков. Длина каждого элементарного отрезка. Точки деления будут: x₀=a; x₁=a+h; x₂=a+2Ч h,., x_n-1=a+ (n-1) Ч h; x_n=b. Эти числа будем называть узлами. Вычислим значения функции f (x) в узлах, обозначим их y₀, y₁, y₂,., y_n. Cталобыть, y₀=f (a), y₁=f (x₁), y₂=f (x₂),., y_n=f (b). Числа y₀, y₁, y₂,., y_n являются ординатами точек графика функции, соответствующих абсциссам x₀, x₁, x₂,., x_n. Площадь криволинейной трапеции приближенно заменяется площадью многоугольника, составленного из n прямоугольников. Таким образом, вычисление определенного интеграла сводится к нахождению суммы n элементарных прямоугольников.

4.3 Формула средних прямоугольников

Рис.2

4.4 Метод правых прямоугольников

Разделим отрезок [a; b] на n равных частей, т. е. на n элементарных отрезков. Длина каждого элементарного отрезка. Точки деления будут: x₀=a; x₁=a+h; x₂=a+2Ч h,., x_n-1=a+ (n-1) Ч h; x_n=b. Эти числа будем называть узлами. Вычислим значения функции f (x) в узлах, обозначим их y₀, y₁, y₂,., y_n. Cталобыть, y₀=f (a), y₁=f (x₁), y₂=f (x₂),., y_n=f (b). Числа y₀, y₁, y₂,., y_n являются ординатами точек графика функции, соответствующих абсциссам x₀, x₁, x₂,., x_n. Площадь криволинейной трапеции приближенно заменяется площадью многоугольника, составленного из n прямоугольников. Таким образом, вычисление определенного интеграла сводится к нахождению суммы n элементарных прямоугольников.

4.5 Формула правых прямоугольников

Рис.3

4.6 Метод Симпсона

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

Разобьем отрезок интегрирования [a; b] на 2Ч n равных частей длины. Обозначим точки разбиения x₀=a; x₁=x₀+h,., x_i=x₀+iЧ h,., x_2n=b. Значения функции f в точках x_i обозначим y_i, т. е. y_i=f (x_i). Тогда согласно методу Симпсона

Рис. 4.

4.7 Метод трапеций

Разделим отрезок [a; b] на n равных частей, т. е. на n элементарных отрезков. Длина каждого элементарного отрезка. Точки деления будут: x₀=a; x₁=a+h; x₂=a+2Ч h,., x_n-1=a+ (n-1) Ч h; x_n=b. Эти числа будем называть узлами. Вычислим значения функции f (x) в узлах, обозначим их y₀, y₁, y₂,., y_n. Cталобыть, y₀=f (a), y₁=f (x₁), y₂=f (x₂),., y_n=f (b). Числа y₀, y₁, y₂,., y_n являются ординатами точек графика функции, соответствующих абсциссам x₀, x₁, x₂,., x_n

Формула трапеций:

Формула означает, что площадь криволинейной трапеции заменяется площадью многоугольника, составленного из n трапеций (рис.5); при этом кривая заменяется вписанной в нее ломаной.

Рис.5

Рис.6

5. Постановка задачи Коши

Определение. Задача нахождения частного решения дифференциального уравнения, удовлетворяющего заданному начальному условию, называется задачей Коши.

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

Определение. Решить задачу Коши для уравнения y'=f (x, y) (6.1) — это значит найти решение уравнения y'=f (x, y) в виде функции у (х), удовлетворяющей начальному условию у (х0) =у0

Геометрически это означает, что требуется найти интегральную кривую у=у (х), проходящую через заданную точку M0 (x0,y0) при выполнении равенства (6.1).

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

Например дифференциальное уравнение у'=у2+х2 не имеет аналитического решения.

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

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

6. Разностные схемы Эйлера

В основе метода Эйлера лежит идея графического построения решения дифференциального уравнения. Однако этот метод дает одновременно и способ нахождения искомой функции в табличной форме.

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

Где, x_i=x₀+iЧ h, — шаг таблицы.

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

y-y₀=f (x₀, y₀) Ч (x-x₀)

y=y₀+f (x₀, y₀) Ч (x-x₀)

если x=x₁, то

y₁=y₀+f (x₀, y₀) Ч (x₁-x₀)

y₁=y₀+hЧ f (x₀, y₀)

D y₀=hЧ f (x₀, y₀)

если x=x₂, то

y₂₌y₁+f (x₁, y₁) Ч (x₂-x₁)

y₂=y₁+hЧ f (x₁, y₁)

D y₁=hЧ f (x₁, y₁)

…

если x=x_i+1, то

y_i+1=y_i+hЧf (x_i, y_i)

Dy_i=hЧ f (x_i, y_i)

Таким образом, получение таблицы значений искомой функции у (х) по методу Эйлера заключается в циклическом применении пары формул:

Dy_k=hЧ f (x_k, y_k)

y_k+1=y_k+Dy_k

где k=0, 1, 2, …, n

Геометрически эти формулы означают, что на отрезке [x_i, x_i+1] интегральная кривая заменяется отрезком касательной к кривой (см. рис. 7, рис.8).

Рис.7 Рис.8

7. Метод Рунге-Кутта второго порядка (Метод Эйлера-Коши)

Рассмотрим метод Рунге-Кутта второго порядка. В этом методе величины y_i+1 вычисляются по следующим формулам:

y_i+1 = y_i + Dy_i

Dy_i=D y_i1+D y_i2

См. рис.1

Тогда .

Обозначим

тогда

Геометрически это означает, что мы определяем направление интегральной кривой в исходной точке (x_i, y_i) и во вспомогательной точке, а в качестве окончательного берем среднее из этих направлений.

8. Метод Рунге-Кутта четвертого порядка

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

В этом методе величины y_i+1 вычисляются по следующим формулам:

y_i+1 = y_i + Dy_i

Dy_i=hЧ (k₁ + 2k₂ + 2k₃ + k₄) /6, i = 0, 1,.

k₁ = f (x_i, y_i),

k₂ = f (x_i+h/2, y_i+hЧ k₁/2),

k₃ = f (x_i+h/2, y_i+hЧ k₂/2),

k₄ = f (x_i+h, y_i+hЧ k₃).

9. Практическая часть

1 Численная реализация решения системы дифференциальных уравнений

1.1 Реализация в программе MathCad. Метода Рунге-Кутта.

Метод Рунге-Кутта

Таблица 1. Результаты работы программы MathCad.

1.2 Реализация в программе MathCad. Метод Эйлера 1-я модификация

Таблица 2. Результаты работы программы MathCad.

Рисунок 2 — График зависимости силы тока от времени

Рисунок 3 — График Зависимости напряжения от времени

1.3 Реализация решения на языке программирования высокого уровня (C++) методом Эйлера (1 модификация).

Блок-схема

Программный код

#include

#include

#include

#include

const double PI = 3.1 415 926;

double f, PhiA, Phi, E0, LA, L, CA, C, t0, t1, t2, h;

double R [7];

int n, i;

// Функция гармонического источника тока, учитывающая переключение ключа

double E (double t)

{

doubleResultE;

if (t < t1)

{

ResultE = E0 + E0 * sin (2 * PI * f* t + Phi);

}

else

{

ResultE = 0;

}

returnResultE;

}

// Производная от тока I по времени t

double FDI (double t, double I, double U)

{

double A, B, D, ResultFDI;

A = R / (R + R [2]);

B = (R * R [2]) / (R + R [2]);

D = (R + R [6]) / (R + R + R [6]);

ResultFDI = (1/L) * ((A * E (t)) — (I * (R + (R * D) + B)) — (U * D));

returnResultFDI;

}

// Производная от напряжения U по времени t

double FDU (double t, double I, double U)

{

double A, B, ResultFDU;

A = (R + R [6]) / (R + R + R [6]);

B = 1/ (R + R + R [6]);

ResultFDU = (1/C) * ((I * A) — (U * B));

returnResultFDU;

}

void main ()

{

FILE *FOUT;

double xt1, xt2, pI, I1, I2, pU, U1, U2;

FOUT = fopen («C: \student\Result08. txt», «a»);

// Вводпараметровцепи

cout<< «********** METOg 3uLEPA (3 mog.) **********n» ;

cout<< «BBEguTE E0 (B): n» ;

cin>> E0;

cout<< «BBEguTE L (MuLurEHPu): n» ;

cin>> LA;

L = LA / 1000;

cout<< «BBEguTE C (MuKPOqpAPAgbI): n» ;

cin>> CA;

C = CA / 1 000 000;

cout<< «BBEguTE t0 (c): n» ;

cin>> t0;

cout<< «BBEguTE t1 (c); n» ;

cin>> t1;

cout<< «BBEguTE t2 (c): n» ;

cin>> t2;

cout<< «BBEguTE n: n» ;

cin>> n;

h = (t2 — t0) / n;

cout<< «War h = «<< h << «n» ;

for (i=1; i<=6; i++)

{

cout<< «BBEguTE R («<< «) (OM): n» ;

cin>> R [i];

}

cout<< «BBEguTE f (Gertz): n» ;

cin>> f;

cout<< «BBEguTE Phi (* PI PAguAH): n» ;

cin>>PhiA;

Phi = PhiA * PI;

// Определение начальных условий

cout<< «BBEguTE I (0), (A): n» ;

cin>> I1;

cout<< «BBEguTE U (0), (B): n» ;

cin>> U1;

xt1 = t0;

cout<< «n» ;

// Запись введенных данных в файл

fprintf (FOUT, «********** Метод Эйлера (3 мод.) **********n»);

fprintf (FOUT, «Параметры цепи: n»);

fprintf (FOUT, «E0 = %g (B) n», E0);

fprintf (FOUT, «L = %g (мГн) n», LA);

fprintf (FOUT, «C = %g (мкФ) n», CA);

fprintf (FOUT, «t0 = %g © nt1 = %g © nt2 = %g © n», t0, t1, t2);

fprintf (FOUT, «n = %dn», n);

fprintf (FOUT, «h = %gn», h);

for (i=1; i<=6; i++)

{

fprintf (FOUT, «R (%d) = %g (Ом) n», i, R [i]);

}

fprintf (FOUT, «f = %g (Гц) n», f);

fprintf (FOUT, «ф = %g*п = %g (рад.) n», PhiA, Phi);

fprintf (FOUT, «nТаблица значений: n»);

fprintf (FOUT, «itttItUn»);

fprintf (FOUT, «0t%gt%gt%gn», xt1, I1, U1);

// Цикл расчета значений функций 3 модификацией метода Эйлера

for (i=1; i<=n; i++)

{

// Вычислениезначений

xt2 = xt1 + h;

pI = I1 + h * FDI (xt1, I1, U1);

pU = U1 + h * FDU (xt1, I1, U1);

I2 = I1 + (h / 2) * (FDI (xt1, I1, U1) + FDI (xt2, pI, pU));

U2 = U1 + (h / 2) * (FDU (xt1, I1, U1) + FDU (xt2, pI, pU));

// Вывод полученных значений на экран

cout<< «i = «<< «tt = «<<< «tI = «<<< «tU = «<<< «n» ;

// Вывод полученных значений в файл

fprintf (FOUT, «%dt%gt%gt%gn», i, xt2, I2, U2);

// Смещение значений для следующей итерации

xt1 = xt2;

I1 = I2;

U1 = U2;

}

cout<< «nnn» ;

fprintf (FOUT, «nРасчетоконченnnn»);

fclose (FOUT);

}

Анализ результатов.

Численное решение системы дифференциальных уравнений в программах MathCad и C++ было реализовано с помощью метода Эйлера 1-й модификации.

Численные результаты двух программ совпадают.

Получены графики зависимости силы тока и напряжения от времени согласно системе уравнений и файл данных, содержащий дискретные зависимости силы тока и напряжения от времени на интервале .

2. Решение задачи аппроксимации зависимости I (t)

2.1 Реализация в пакете Excel

Таблица 4. Результаты работы программы Excel

Рисунок 7 — график зависимости I (t) и аппроксимирующая функция на первом участке Рисунок 8 — график зависимости I (t) и аппроксимирующая функция на втором участке Рисунок 9 — график зависимости I (t) и аппроксимирующая функция на третьем участке Рисунок 10 — график зависимости I (t) и аппроксимирующая функция на трех участках Анализ результатов Решение задачи аппроксимации было проверено в MathCAD и в MicrosoftExcel 2003. Точки для аппроксимации брались из пункта (MathCAD). В среде MathCAD был реализован метод наименьших квадратов, при этом участок был разбит на три промежутка для получения более точной аппроксимирующей функции. В Excel аппроксимация проведена с помощью линий тренда.

2.2 Реализация в программе MathCad. Метод наименьших квадратов

Первый участок от t1=0 до t2=0.0002

Второй участок от t=0.0002 до t=0.0006

программный код численное интегрирование

Третий участок от t3=0.0006 до t4=0.0099

Совмещение полученных аппроксимирующих функций

Анализ результатов.

Решение задачи аппроксимации было проведено в программе MathCad методом наименьших квадратов и в пакете Excel с использованием мастера диаграмм с выводом уравнения линии тренда. Результаты двух программ совпадают. Получена аналитическая формула зависимости силы тока от времени на интервале 0

3 Численное интегрирование

3.1 Реализация в программе С++. Метод левых прямоугольников.

3.1.1 Блок-схема Рисунок 16 — Блок-схема. Метод левых прямоугольников

3.1.2 Код программы

#include" stdafx. h"

#include

#include

#include

usingnamespacestd;

// Квадрат аналитической функции I (t), полученной в пакете Excel

double I (double t)

{doubleIt, Res;

if (t<0.0004)

{It=4* (1E+9) * (pow (t, 3)) — 5* (1E+6) *pow (t, 2) +1873.5*t+0.0002; }

else

if (t<0.001)

{It =2* (1E+8) *pow (t, 3) — 513 156*pow (t, 2) — 346.61*t+0.1905; }

else

{It =-2* (1E+6) *pow (t, 3) +19 110*pow (t, 2) — 99.733*t+0.3122; }

Res = pow (It, 2);

return Res; }

void main ()

{double R4, T1,T2,h, S, Integ, Q;

inti, n;

cout<<<" metodLevihpryamougolnikov" <

cout<<" Vveditesoprotivltnie R4 (OM)" <

cin>>R4;

cout<<" Vvedite T1 (c): «<

cin>>T1;

cout<<" Vvedite T2 (c): «<

cin>>T2;

cout<<" Vveditekolichectvorasbienii: «<

cin>> n;

h= (T2 — T1) /n;

cout<<" h = «<

S=0;

for (i=1; i<= (n-1); i++)

{S=S+I (T1+h*i); }

Integ = h*S;

Q = R4*Integ;

cout<<" Integ = «<

cout<<" Kolichestvoteploti Q ="<

Метод правых прямоугольников.

#include" stdafx. h"

#include

#include

#include

usingnamespacestd;

// Квадрат аналитической функции I (t), полученной в пакете Excel

double I (double t)

{double It, Res;

if (t<0.0004)

{It=4* (1E+9) * (pow (t, 3)) — 5* (1E+6) *pow (t, 2) +1873.5*t+0.0002; }

else

if (t<0.001)

{It =2* (1E+8) *pow (t, 3) — 513 156*pow (t, 2) — 346.61*t+0.1905; }

else

{It =-2* (1E+6) *pow (t, 3) +19 110*pow (t, 2) — 99.733*t+0.3122; }

Res = pow (It, 2);

return Res; }

void main ()

{double R4, T1,T2,h, S, Integ, Q;

int i, n;

cout<<<" metod Pravih pryamougolnikov" <

cout<<" Vvedite soprotivltnie R4 (OM)" <

cin>>R4;

cout<<" Vvedite T1 (c): «<

cin>>T1;

cout<<" Vvedite T2 (c): «<

cin>>T2;

cout<<" Vvedite kolichectvo rasbienii: «<

cin>> n;

h= (T2 — T1) /n;

cout<<" h = «<

S=0;

for (i=1; i<=n; i++)

{

S=S+I (T1+h*i); }

Integ = h*S;

Q = R4*Integ;

cout<<" Integ = «<

cout <<" Kolichestvo teploti Q ="<

Анализ результатов.

Численное решение системы дифференциальных уравнений в программах MathCad и C++ было реализовано с помощью метода Эйлера 1-й модификации. Численные результаты двух программ совпадают. Получены графики зависимости силы тока и напряжения от времени согласно системе уравнений и файл данных, содержащий дискретные зависимости силы тока и напряжения от времени на интервале 0.004

Реализация в программе MathCad.

Расчет количества теплоты, выделившейся на резисторе R4

1 Метод трапеций

2 Метод левых прямоугольников:

3 Метод правых прямоугольников:

4 Метод Симпсона:

5 Метод центральных прямоугольников

Вычисление ошибок

Анализ результатов.

Численное интегрирование было реализовано в программе С++ методом центральных прямоугольников и в программе MathCad разными методами численного интегрирования. Результаты двух программ совпадают. Разные методы дают разную точность вычисления, что видно по рассчитанным ошибкам. Наибольшую точность имеет решение, полученное методом левых и правых прямоугольников.

Выводы

В данной работе были проанализированы переходные процессы в электрической цепи переменного тока. Все расчеты проведены с помощью численных методов решения математических задач. Была получена система дифференциальных уравнений для и, которая была решена модифицированным методом Эйлера в программах MathCad и С++. Была решена задача аппроксимации полученной дискретной зависимости при помощи пакета Excel и MathCad, в результате было получено аналитическое уравнение зависимости. Используя это уравнение, при помощи численных методов интегрирования, наиболее точным из которых оказался метод левых и правых прямоугольников, было найдено количество теплоты, выделяемое на резисторе .

Задание 1. Сравнение результатов

Таблица 6. MathCAD. Метод Эйлера1-ая модификация

Таблица 7. C++ метода Эйлера2-й модификации

Таблица 8. MathCAD. Метод Рунге-Кутта

Задание 2. Сравнение результатов Таблица 9. Значения теплоты, полученные в MathCAD

1. Численное моделирование и анализ переходных процессов в электрической цепи. Метод. разработка по выполнению курсовой работы по информатике для студентов технических специальностей дневной формы обучения/НГТУ; Сост. С. Н. Митяков, Т. В. Моругина, М. Н. Потапова, Т. А. Факеева. Н. Новгород, 2004. — 12 с.

2. Численные методы анализа. Приближение функции, дифференциальные и интегральные уравнения / Б. П. Демидович, И. А. Марон, Э. З. Шувалова. — СПб.: Издательство «Лань», 2005. — 655 с.

3. Самоучитель MathCAD 11. /Кирьянов Д.В. — СПб.: БХВ-Петербург, 2003. — 560 с. ил.

4. C++. Программирование на языке высокого уровня. Учебник для вузов/Павловская Т.А. — Питер, 2004. — 393 с.

5. Информатика и информационные технологии. Учебное пособие/И.Г. Лесничая, И. В. Миссинг.2-е изд. — М.: Изд-воЭксмо, 2008. — 544 с. (Высшее экономическое образование)