Задачи, алгоритмы, программы (стр. 3 )

Рис. 3.8. Конвекция двумерного газа.

Задача 3.9. Рост некоторых природных объектов (кораллов, снежных хлопьев) происходит в результате случайного присоединения частей. Этот процесс называется агрегацией с ограничением диффузии. Возникает фрактал, имеющий древовидную структуру и называющийся дендритом. Напишите программу, моделирующую рост дендрита.

Промоделируем хаотическое движение молекул некоторого вещества в растворе. Все молекулы находятся в состоянии 0. Когда какая-нибудь молекула достигнет ячейки с заравочной частицей, она перейдет в состояние 1, остановится и изменит свой цвет. Если другая частица приблизится к молекуле, находящейся в состоянии 1, она как бы осядет, присоединится к ней, тоже перейдя в состояние 1. Частицы, находящиеся в состоянии 0, продолжают хаотически двигаться и оседать на частицах, находящихся в состоянии 1. В результате растет фрактал, имеющий древовидную структуру. Частицы, обуславливающие его рост, не могут проникнуть в глубь образующегося дендрита, они оседают на его периферии. Используется программа ПР - 3.9, результат ее работы представлен на рис. 3.9.2. В ней используются периодические граничные условия: когда частица выходит через левую границу квадратной ячейки, она входит в нее через правую границу; когда частица выходит через нижнюю границу квадратной ячейки, она входит в нее через верхнюю границу и т. д.

{$N+} uses dos, crt, graph; { ПР - 3.9 } const N=400; var Fx, Fy, x,y, vx, vy: array[1..N] of Longint; h: array[1..N] of boolean; r, t,Gd, Gm, i,j: integer; k, ax, ay, F,l: single; hh: string; v: real; label m0; Procedure Raschet; label m1; begin For i:=1 to N do For j:=1 to N do If h[i]=false then begin l:=sqrt(sqr(x[i]-x[j])+sqr(y[i]-y[j])); If (h[j]=true)and(l<6) then h[i]:=true; end; end; Procedure Nach_uslov; begin For i:=1 to N do begin x[i]:=100+round(random(300)); y[i]:=100+round(random(300)); end; end; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); Nach_uslov; Randomize; setbkcolor(15); Repeat Raschet; For i:=1 to N do begin If h[i]=true then goto m0; x[i]:=x[i]+round(random(600)/100)-3; y[i]:=y[i]+round(random(600)/100)-3; If x[i]<100 then x[i]:=400; If x[i]>400 then x[i]:=100; If y[i]>400 then y[i]:=100; If y[i]<100 then y[i]:=400; If (abs(x[i]-250)<4)and(abs(y[i]-250)<4)then h[i]:=true; m0: end; If t mod 20=0 then begin cleardevice; For i:=1 to N do begin If h[i]=true then setcolor(blue) else setcolor(red); circle(round(x[i]),round(y[i]),2); circle(round(x[i]),round(y[i]),1); end; end; inc(t); until KeyPressed; Repeat until keypressed; CloseGraph; END.

Рис. 3.9.1. Рост дендрита.

Рис. 3.9.2. Дендрит имеет фрактальную структуру.

4. МЕХАНИКА ТВЕРДОГО ТЕЛА

Задача 4.1. Стержень устанавливают на горизонтальную поверхность так, чтобы он образовывал некоторый угол с вертикалью и отпускают. Необходимо рассчитать координаты его концов при падении в последовательные моменты времени.

Рис. 4.1.1. Падение стержня на горизонтальную поверхность.

Рассмотрим систему, состоящую из двух материальных точек m1 и m2, связанных между собой невесомым упругим стержнем. Длина стержня в недеформированном состоянии равна l0, при его сжатии возникают силы упругости F1 и F2. Будем считать, что нижний конец стержня A скользит по горизонтальной поверхности, не отрываясь от нее (y1=0). При этом на него действует сила вязкого трения FТР, пропорциональная скорости и направленная в сторону противоположную движению. Проекции сил, действующих на материальные точки m1 и m2 вычисляются из формул: В программе ПР - 4.1 в последовательные моменты времени пересчитывают силы, действующие на материальные точки, их ускорения, скорости и координаты. Результаты моделирования падения стержня при различных значениях коэффициента сопротивления r представлены на рис. 4.1.2. При необходимости можно рассчитать траекторию движения центра масс стержня, зависимость его угла наклона от времени и т. д.

uses dos, crt, graph; { ПР - 4.1 } const m1=1; m2=2; k=100; r=0.02; l0=300; dt=0.001; var l, x1,x2,y1,y2,vx1,vy1,vx2,vy2,F, Fx1,Fy1,Fx2,Fy2 : real; Gd, Gm, n : integer; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); x1:=200; y1:=0; x2:=201; y2:=300; line(0,452,640,452); Repeat l:=sqrt(sqr(x1-x2)+sqr(y1-y2)); F:=k*(l0-l); Fx1:=-F*(x2-x1)/l-r*vx1; Fx2:=F*(x2-x1)/l; Fy2:=F*(y2-y1)/l-9.8*m2; vx1:=vx1+Fx1/m1*dt; vy1:=vy1+Fy1/m1*dt; vx2:=vx2+Fx2/m2*dt; vy2:=vy2+Fy2/m2*dt; x1:=x1+vx1*dt; x2:=x2+vx2*dt; y2:=y2+vy2*dt; if y2<1 then begin vx2:=0; vy2:=0; end; if n mod 500=0 then begin line(20+round(x1),450-round(y1),20+round(x2),450-round(y2)); line(21+round(x1),451-round(y1),21+round(x2),451-round(y2)); circle(20+round(x2),450-round(y2),3); end; n:=n+1; until KeyPressed; END.

Рис. 4.1.2. Результаты моделирования падения стержня.

Задача 4.2. Лестницу прислоняют к стене и отпускают. Напишите программу, моделирующую падение лестницы, если известно, что ее концы с трением скользят по поверхностям стены и пола, не отрываясь от них.

Рис. 4.2. Падение стержня на горизонтальную поверхность.

Вместо лестницы рассмотрим систему двух материальных точек массами m1 и m2, соединенных невесомым упругим стержнем жесткостью k и длиной l0. Эти точки как бы скользят по вертикальной и горизонтальной направляющим, при этом на них действует сила вязкого трения, прямо пропорциональная скорости и направленная в противоположную сторону. Найдем действующие на точки силы: Программа ПР - 4.2 содержит цикл по времени в котором вычисляются силы, действующие на материальные точки, определяются их ускорения, скорости и координаты, а также координаты центра масс. На экране монитора изображается положение лестницы в последовательные моменты времени.

uses dos, crt, graph; { ПР - 4.2 } const m1=3; m2=1; k=1000; r=0.1; dt=0.0004; var x1,x2,y1,y2,xc, yc, vx1,vy1,vx2,vy2,l, l0,F, Fx1,Fy1,Fx2,Fy2 : real; Gd, Gm, n : integer; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); x1:=100; y1:=300; x2:=150; y2:=0; l0:=sqrt(sqr(x1-x2)+sqr(y1-y2)); rectangle(100,450,500,100); Repeat l:=sqrt(sqr(x1-x2)+sqr(y1-y2)); F:=k*(l0-l); Fy1:=F*(y1-y2)/l-9.8*m1-r*vy1; Fx2:=F*(x2-x1)/l-r*vx2; vy1:=vy1+Fy1/m1*dt; vx2:=vx2+Fx2/m2*dt; x2:=x2+vx2*dt; y1:=y1+vy1*dt; xc:=(m1*x1+m2*x2)/(m1+m2); yc:=m1*y1/(m1+m2); if y1<1 then begin vx1:=0; vy1:=0; end; delay(10); if n mod 2000=0 then begin setcolor(15); line(round(x1),450-round(y1),round(x2),450-round(y2)); line(round(x1+1),451-round(y1),round(x2+1),451-round(y2)); setcolor(12); circle(round(x1),450-round(y1),2); circle(round(x2),450-round(y2),2); circle(round(xc),450-round(yc),3); circle(round(xc),450-round(yc),4); end; n:=n+1; until KeyPressed; END.

Рис. 4.2. Падение стержня на горизонтальную поверхность.

Задача 4.3. Промоделируйте плоское движение тела в поле тяжести, при котором его центр масс движется по кривой, а само тело вращается вокруг центра масс. Используйте модель состоящую из двух материальных точек, связанных упругим стержнем.

Рис. 4.3.1. К вычислению массы, объема и координат центра масс тела.

Рассмотрим систему из двух материальных точек, соединенных упругим стержнем. Проекции сил, действующих на точки, могут быть вычислены по формулам: <="" p=""> Программа ПР - 4.3 содержит цикл по времени, в котором вычисляются действующие силы, ускорения, скорости и координаты точек. На экране строятся положения стержня в последовательные моменты времени. Программа позволяет промоделировать: 1) отскок стержня от стены и пола; 2) движение системы в случае, когда массы точек сильно отличаются; 3) движение системы при не очень большой жесткости стержня.

uses dos, crt, graph; { ПР - 4.3 } const m1=0.2; m2=0.1; k=1000; r=0.005; b=0.6; l0=40; dt=0.001; var l, x1,x2,y1,y2,vx1,vy1,vx2,vy2,F, Fx1,Fy1,Fx2,Fy2 : real; Gd, Gm, n : integer; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); x1:=0; y1:=300; x2:=0; y2:=260; vx1:=50; vy1:=20; rectangle(2,2,572,442); Repeat l:=sqrt(sqr(x1-x2)+sqr(y1-y2)); F:=k*(l-l0); Fx1:=F*(x2-x1)/l-r*vx1; Fy1:=F*(y2-y1)/l-r*vy1-0.98*m1; Fx2:=-F*(x2-x1)/l-r*vx2; Fy2:=-F*(y2-y1)/l-r*vy2-0.98*m2; vx1:=vx1+Fx1/m1*dt; vy1:=vy1+Fy1/m1*dt; vx2:=vx2+Fx2/m2*dt; vy2:=vy2+Fy2/m2*dt; x1:=x1+vx1*dt; y1:=y1+vy1*dt; x2:=x2+vx2*dt; y2:=y2+vy2*dt; If x1>550 then vx1:=-b*vx1; If x2>550 then vx2:=-b*vx2; If y1<-0 then vy1:=-b*vy1; If y2<-0 then vy2:=-b*vy2; circle(20+round(x2),440-round(y2),1); If n mod 1000=0 then begin line(20+round(x1),440-round(y1),20+round(x2),440-round(y2)); circle(20+round(x1),440-round(y1),3); end; inc(n); until KeyPressed; END.

Рис. 4.3.2. Движение палки в поле тяжести.

Задача 4.4. На тележке массой m1 подвешен маятник, состоящий из тела массой m2 и нити длиной l. Маятник выводят из положения равновесия и отпускают. В подшипниках тележки действует сила вязкого трения. Напишите программу, моделирующую затухающие колебания системы.

Рис. 4.4.1. Движение тележки с маятником.

Заменим систему "маятник-тележка" системой, состоящей из двух материальных точек m1 и m2, связанных упругим стержнем жесткостью k и длиной l0. Материальная точка m1 способна скользить по горизонтальной линии так, что ее координата y2остается постоянной. При этом на нее действует сила вязкого трения, направленная противоположно скорости и пропорциональная ее величине. Проекции сил, действующих на точки системы, вычисляются по формулам: В предлагаемой программе ПР - 4.4 осуществляется расчет действующих на точки сил, а также проекций их ускорений, скоростей и координат в последовательные моменты времени. Результаты моделирования представлены на рис. 4.4.

uses dos, crt, graph; { ПР - 4.4 } const m1=0.1; m2=0.2; k=100; r=0.01; dt=0.001; var l0,l, x1,x2,y1,y2,vx1,vy1,vx2,vy2,F, Fx1,Fy1,Fx2,Fy2 : real; Gd, Gm : integer; n: longint; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); x1:=0; y1:=0; x2:=300; y2:=-150; l0:=sqrt(sqr(x1-x2)+sqr(y1-y2)); Repeat l:=sqrt(sqr(x1-x2)+sqr(y1-y2)); F:=k*(l-l0); Fx1:=F*(x2-x1)/l-r*vx1; Fx2:=-F*(x2-x1)/l; Fy2:=-F*(y2-y1)/l-0.5*m2; vx1:=vx1+Fx1/m1*dt; vx2:=vx2+Fx2/m2*dt; vy2:=vy2+Fy2/m2*dt; x1:=x1+vx1*dt; x2:=x2+vx2*dt; y2:=y2+vy2*dt; circle(250+round(x2),50-round(y2),1); if n mod 8000=0 then begin line(250+round(x1),50-round(y1),250+round(x2),50-round(y2)); line(251+round(x1),51-round(y1),251+round(x2),51-round(y2)); circle(250+round(x1),50-round(y1),2); circle(250+round(x2),50-round(y2),2); circle(250+round(x2),50-round(y2),3); end; inc(n); until KeyPressed; END.

Рис. 4.4.2. Колебания маятника на тележке.

Задача 4.5. На горизонтальной поверхности покоится кольцо (труба), к внутренней стороне которого прикреплен груз. Расстояние от оси кольца до его центра масс известно. Кольцо смещают из положения равновесия от отпускают. Изучите: 1) колебания кольца относительно положения равновесия; 2) движение кольца после того, как ему сообщили начальную скорость.

Рис. 4.5.1. Колебания маятника на тележке.

Необходимо рассчитать расстояние от центра кольца O до центра масс C, момент M силы тяжести, момент инерции I относительно мгновенной оси вращения A. Для этого используют формулы: Зная момент силы и момент инерции определяют угловое ускорение тела в последовательные моменты времени, вычисляют угловую скорость и угол поворота, горизонтальную координату центра кольца. При этом используется программа ПР-1, результаты представлены на рис. 4.5.2, 4.5.3. Чтобы промоделировать качение кольца со смещенным центром масс, следует задать начальную скорость (рис. 4.5.4).

uses dos, crt, graph; { ПР - 4.5 } const dt=0.005; g=9.8; R=60; l=55; var x, w,fi, mt, mgr, lc, M,I, dd : real; Gd, Gm : integer; k: longint; BEGIN Gd:=Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); mt:=1; mgr:=10; x:=-10; fi:=3.0;w:=0;{-0.1;} Repeat inc(k); lc:=mgr*l/(mgr+mt); M:=(mgr+mt)*g*lc*sin(fi); dd:=R*R+l*l-2*R*l*cos(fi); I:=mgr*dd+2*mt*R*R; w:=w-M/I*dt; fi:=fi+w*dt; x:=x-R*w*dt; { circle(10+round(k/100),300-round(0.4*x),1);} circle(100+round(x+R*sin(fi)),300+round(R*cos(fi)),1); if k mod 500=0 then begin {cleardevice;} line(0,300+R,640,300+R); line(0,301+R,640,301+R); line(100+round(x),300,100+round(x+R*sin(fi)),300+round(R*cos(fi))); circle(100+round(x),300,R); circle(100+round(x),300,R-1); end; until (k>7000)or(KeyPressed); Repeat until KeyPressed; CloseGraph; END.

Рис. 4.5.2. Зависимость угла поворота от времени.

Рис. 4.5.3. Колебания кольца на горизонтальной поверхности.

Рис. 4.5.4. Качение кольца со смещенным центром масс.

Задача 4.6. Брусок в форме прямоугольного параллелепипеда плавает на поверхности жидкости. Координаты центра масс, плотности жидкости и бруска известны. Определите расположение бруска относительно поверхности жидкости.

Рис. 4.6.1. К расчету положения тела.

Расположение бруска относительно поверхности воды однозначно определяется величинами d и φ (рис. 4.6.1). Пусть a=3, b=1. Задача решается так. Положим, что φ=0. Будем увеличивать d от -3 с некоторым шагом, каждый раз вычисляя объем погруженной части бруска и определяя силу Архимеда. Плотность бруска ρ меньше плотности жидкости ρ0, поэтому сила Архимеда превысит силу тяжести, когда тело не полностью погрузится в жидкость. Нарисуем поверхность жидкости AB и определим положение центра плавучести P. Брусок будет находиться в положении устойчивого равновесия тогда, когда центр масс C находится точно под центром плавучести P, то есть углы φ и φ1 равны. Если это условие не выполняется, будем увеличивать угол φ с некоторым шагом, каждый раз определяя величину d, положение центра плавучести P и проверяя равенство углов φ и φ1 (оно может выполняться с небольшой погрешностью). Когда углы φ и φ1 окажутся равными, программы должна вывести результат вычислений. Для нахождения объема погруженной части тела Vp необходимо найти площадь S фигуры, выделенной точками. Для этого используется метод прямоугольников: For i:=-100 to 100 do begin x:=dx*i; y:=sin(fi)/cos(fi)*x+d; If y>b then y:=b; If y<-b then y:=-b; S:=S+(y+b)*dx; end; Чтобы найти координаты xpl и ypl центра плавучести, внутрь прямоугольника случайным образом бросим 100000 точек и подсчитаем количество точек, попавших в заштрихованную область S, соответсвующую погруженной части тела. Используется следующий фрагмент программы: Repeat x:=(random(200)/100-1)*a; y:=(random(200)/100-1)*b; yl:=sin(fi)/cos(fi)*x+d; inc(N); if y100000; xpl:=sx/(nn+0.01); ypl:=sy/(nn+0.01); Результаты работы программы ПР - 4.6 представлены на рис. 4.6.2 и 4.6.3. Программа рисует прямоугольное сечение бруска и определяет положение поверхности жидкости AB и центра плавучести. Сила Архимеда направлена вдоль вертикали MN. В состоянии устойчивого равновесия центр масс C и центр плавучести P лежат на одной вертикали MN. Изображение, получающееся на экране монитора, представлено в левой части рис. 4.6.2 и 4.6.3.

{$N+} uses dos, crt, graph; const rho=900; rho0=1000; g=9.8; { ПР - 4.6 } Ms=50; a=3; b=1; c=1; xc=1.2; yc=-0.5; pi=3.1415; var x, y,dx, d,yl, Vp, V,S, sx, sy, xpl, ypl, Fa, Ft, m,fi, fi1: single; i, j,k, Gd, Gm: integer; n, nn: Longint; BEGIN Gd:= Detect; InitGraph(Gd, Gm, 'c:\bp\bgi'); Randomize; dx:=a/100; Repeat inc(j); fi:=pi/300*j; d:=-3; Repeat S:=0; {cleardevice;} For i:=-100 to 100 do begin x:=dx*i; y:=sin(fi)/cos(fi)*x+d; If y>b then y:=b; If y<-b then y:=-b; S:=S+(y+b)*dx; end; V:=4*a*b*c; Vp:=S*c; Fa:=rho0*Vp*g; Ft:=rho*4*a*b*c*g; d:=d+0.01; until Fa>Ft; rectangle(320-Ms*a,240-Ms*b,320+Ms*a,240+Ms*b); rectangle(318-Ms*a,238-Ms*b,322+Ms*a,242+Ms*b); circle(round(320+Ms*xc),round(240-Ms*yc),5); N:=0; nn:=0; sx:=0; sy:=0; Repeat x:=(random(200)/100-1)*a; y:=(random(200)/100-1)*b; yl:=sin(fi)/cos(fi)*x+d; inc(N); If y<yl then begin inc(nn); sx:=sx+x; sy:=sy+y; end; until N>100000; xpl:=sx/(nn+0.01); ypl:=sy/(nn+0.01); fi1:=arctan((xc-xpl)/abs(yc-ypl)); circle(round(320+Ms*xpl),round(240-Ms*ypl),2); circle(round(320+Ms*xpl),round(240-Ms*ypl),3); line(round(320+Ms*(xpl-5*sin(fi))),round(240- Ms*(ypl+5*cos(fi))), round(320+Ms*(xpl+5*sin(fi))), round(240-Ms*(ypl-5*cos(fi)))); line(320+Ms*5,240-round(Ms*sin(fi)/cos(fi)*5+Ms*d), 320-Ms*5,240-round(-Ms*sin(fi)/cos(fi)*5+Ms*d)); until fi1<fi; Repeat until KeyPressed; CloseGraph; END.

Рис. 4.6.2. Устойчивое состояние плавающего тела.

Рис. 4.6.3. Устойчивое состояние плавающего тела.

5. РАСЧЕТ ЭЛЕКТРИЧЕСКИХ ЦЕПЕЙ

Задача 5.1. Электрическая цепь (рис. 5.1) состоит из источника постоянной ЭДС E=200 В, резистора R1=400 Ом и нелинейного элемента, вольт--амперная характеристика которого может быть апроксимирована многочленом третьей степени i(u)=2+5u-0,18u2+0,002u3. Рассчитайте ток в цепи.

Рис. 5.1. Расчет нелинейной цепи.

Вольт--амперные характеристики обоих элементов изображены на рис. 5.1 (ВАХ резистора повернута вокруг оси ординат и смещена на E). Напряжение на нелинейном элементе равно U=E-IR1, поэтому для нахождения искомого тока I следует решить нелинейное алгебраическое уравнение: Чтобы решить уравнение необходимо определить количество корней и локализовать их, то есть установить интервалы, внутри которых находится каждый корень. Воспользуемся графическим методом, для этого построим график функции и найдем точку его пересечения с осью абсцисс, для которой F(I)=0. Она лежит в интервале 0,2 -- 0,4 мА. Самым простым способом решения подобных уравнений является табулирование функции F(I). Это может быть реализовано с помощью алгоритма A-1 (программа ПР - 5.1.1). При запуске программы на экран выводятся те значения тока, для которых функция F(I) по модулю меньше чем ε =0,01, то есть близка к нулю. Это позволяет уменьшить неопределенность в нахождении корня до 0,01, корень уравнения лежит в интервале от 2,960 до 2,965. Ответ: ток примерно равен 2,962 А.

uses crt, graph; const N=200; { ПР - 5.1.1 } var EC, D,DrV, MV, i,j, q1,q2,q3,x1,y1,x2,y2,x3,y3,k: integer; u, tok, tok1,a0,a1,a2,a3,r1,E: real; BEGIN clrscr; a0:=2; a1:=5; a2:=-0.18; a3:=0.002; R1:=400; E:=200; For i:=0 to 2000 do begin tok:=0.0005*i; u:=E-Tok*R1; Tok1:=0.001*(a0+a1*u+a2*u*u+a3*u*u*u); If abs(tok-tok1)<0.01 then writeln(tok, ' ', tok1, ' ', tok-tok1); end; readkey; END.

Другой способ решения уравнения F(x)=0 состоит в использовании метода половинного деления. Задается интервал [a0, b0], содержащий корень. График функции y=F(x) должен пересекать ось абсцисс внутри интервала, на его границах функция имеет противоположны знаки, в чем можно убедиться, проверив условие F(a0)F(b0)<0. Отрезок [a0, b0] делится точкой с=(a+b)/2 попалам. Если F(c)=0, то значение x=c и есть корень уравнения. В противном случае из двух отрезков [a0,c] и [c, b0] выбирается тот, внутри которого график y=F(x) пересекает ось абсцисс. Для этого снова проверяется условие F(c)F(b0)<0, если оно выполняется, то a1=c, b1=b0. После этого осуществляется следующая итерация: отрезок [a1, b1] снова делится попалам и т. д. Так продолжается до тех пор, пока длина отрезка [ai, bi] не окажется меньше заданной точности ε. Используется программа ПР - 5.1.2. Понятно, что ток в цепи не может превышать значение E/R1, поэтому в качестве нулевого приближения к корню выберем интервал [0, E/R1]. В результате последовательности итераций при ε=0,0000001 получается корень уравнения 2,96296. Округляя до четырех значащих цифр, получаем: ток равен 2,963 А.

uses crt; { ПР - 5.1.2 } var a0,a1,a2,a3,r1,E, u,tok, tok1,y1,y2,y3,a, b,c, eps: real; Label metka; Function Funk(Tok: real): real; begin u:=E-Tok*R1; Funk:=0.001*(a0+a1*u+a2*u*u+a3*u*u*u)-Tok; end; BEGIN clrscr; a0:=2; a1:=5; a2:=-0.18; a3:=0.002; R1:=400; E:=200; a:=0; b:=E/R1; eps:=0.0000001; metka: y1:=Funk(a); y2:=Funk(b); c:=(a+b)/2; y3:=Funk(c); writeln(a, b); if (y1*y3)>0 then a:=c else b:=c; If (b-a)>eps then goto metka; writeln('Корень лежит в интервале ', a, b); writeln('0= ',Funk(a)); readkey; END.

Задача 5.2. Напишите программу, рассчитывающую сопротивление цепи, состоящей из бесконечной цепочки резисторов R1=80 Ом и R2=100 Ом (рис. 5.2).

Рис. 5.2. Схема из бесконечной цепочки резисторов.

Найдем сопротивление участка цепи MABN из двух ветвей, участка цепи MACDBN из четырех ветвей, участка цепи MACEFDBN из шести ветвей. Обобщим получающуюся формулу: Программа ПР - 5.2 последовательно вычисляет сопротивление Z1, Z2, Z3, ... и результат выводит на экран. Получающиеся значения стремятся к величине 138 Ом.

uses crt; { ПР - 5.2 } const r1=80; r2=100; var i: integer; z: array [1..200] of real; BEGIN clrscr; i:=1; z[1]:=r1+r2; Repeat writeln(i,' ',z[i]); inc(i); z[i]:=r1+r2*z[i-1]/(r2+z[i-1]); until (i>20)or(KeyPressed); readkey; END.

Задача 5.3. Напишите программу, решающую систему из N линейных алгебраических уравнений.

Рассмотрим алгоритм, позволяющий решать систему из n линейных алгебраических уравнений с n неизвестными, имеющую следующий вид: От этой системы перейдем к матрице из элементов ai, j, где i=1, 2, ... , n и j=1, 2, ... , n+1. Приведем данную матрицу к треугольному виду. Чтобы исключить деление на ноль сложим левые и правые части первого и второго, первого и третьего, второго и четвертого уравнений и т. д. Первый шаг состоит в следующем: разделим все элементы каждой строки матрицы (или каждого уравнения системы) на первый коэффициент ai1, где i=1, 2, ... , n. Тогда все элементы ai1 из первого столбца будут равны 1. Вычтем из второго, третьего и последующих уравнений первое. Матрица коэффициентов приобретет вид: где a'i, j=ai, j/ai1- a1j/a11, где i=1, 2, ... , n и j=1, 2, ... , n+1. После этого сделаем второй шаг: разделим вторую, третью и последующие строки матрица на соответствующие коэффиценты $ai2, после чего вычтем из третьей, четвертой и последующих строк вторую строку. Получим: где a''ij=a'ij/a'i1- a'1j/a'11, где i=1, 2, ... , n и j=1, 2, ... , n+1. После n-ого шага получим матрицу треугольного вида: Искомые значения находятся по формулам: Желательно, чтобы программа осуществляла проверку решения, подставляя найденные значения x1, x2, ... в исходную систему уравнений.

uses crt; const n=6; var i, ii, j,k : integer; { ПР - 5.3 } aa, a: array[1..N, 1..N+1] of real; b, x: array[1..N] of real; BEGIN clrscr; a[1,1]:=13; a[1,2]:=21; a[1,3]:=0; a[1,4]:=-8; a[1,5]:=0; a[1,6]:=0; a[1,7]:=-24; a[2,1]:=0; a[2,2]:=0; a[2,3]:=0; a[2,4]:=8; a[2,5]:=-17;a[2,6]:=11; a[2,7]:=-14; a[3,1]:=0; a[3,2]:=-21;a[3,3]:=15; a[3,4]:=0; a[3,5]:=0; a[3,6]:=-11;a[3,7]:=30; a[4,1]:=1; a[4,2]:=0; a[4,3]:=0; a[4,4]:=1; a[4,5]:=1; a[4,6]:=0; a[4,7]:=0; a[5,1]:=-1; a[5,2]:=1; a[5,3]:=1; a[5,4]:=0; a[5,5]:=0; a[5,6]:=0; a[5,7]:=0; a[6,1]:=0; a[6,2]:=0; a[6,3]:=-1; a[6,4]:=0; a[6,5]:=-1; a[6,6]:=-1; a[6,7]:=0; For i:=1 to N do For j:=1 to N+1 do aa[i, j]:=a[i, j]; For i:=1 to N do For j:=1 to N+1 do a[i, j]:=a[i, j]+a[1,j]+a[N, j]; For i:=2 to N do For j:=1 to N+1 do a[i, j]:=a[i, j]+a[2,j]+a[N-1,j]; For i:=3 to N do For j:=1 to N+1 do a[i, j]:=a[i, j]+a[3,j]+a[4,j]; For i:=1 to N do begin For j:=1 to N+1 do write(a[i, j]:2:3, ' '); writeln; end; writeln; For k:=1 to N do begin For i:=k to N do For j:=1 to N+1 do a[i, N+2-j]:=a[i, N+2-j]/a[i, k]; For i:=k+1 to N do For j:=1 to N+1 do a[i, j]:=a[i, j]-a[k, j]; end; For i:=1 to N do begin For j:=1 to N+1 do write(a[i, j]:2:3, ' '); writeln; end; writeln; delay(13000); For ii:=1 to N do begin i:=N-ii+1; x[i]:=-a[i, N+1]; For j:=i+1 to N do x[i]:=x[i]-a[i, j]*x[j]; end; writeln(' ОТВЕТ '); For i:=1 to N do writeln('x[',i,']= ',x[i]); writeln(' ПРОВЕРКА '); For i:=1 to N do begin For j:=1 to N do b[i]:=b[i]+aa[i, j]*x[j]; writeln(b[i]+aa[i, N+1]:2:7); end; ReadKey; END.

Задача 5.4. Рассчитате цепь постоянного тока, состоящую из нескольких контуров (рис. 5.4). Сопротивления резисторов R1=13 Ом, R2=21 Ом, R3=15 Ом, R4=8 Ом, R5=17 Ом, R6=11 Ом. ЭДС источников E1=5 В, E2=13 В, E3=9 В, E4=-6 В, E5=12 В, E6=-8 В.