Гурский Д., Турбина Е. - Вычисления в MathCad 12 (1077322), страница 89
Текст из файла (страница 89)
Уравнение Лапласатакже применяется для описания полей со сферической симметрией.390 *Глава 11. Вычисление производныхПример 11.14. Найти лапласиан функции и(р, 0, <|))=p2-sin(0)-cos(<|))в сферической системе координатНачнем мы с того, что зададим саму функцию и:и := р • s i n ( 8 ) - c o s ^ )Далее введем выражения, связывающие декартовы координаты х, у, z со сферическими координатами р, 6, ф:X:=p-sin(e)-cos(<|>)Y:«=p-sm(e)-sin(<t>)Z:=pcos(0)Для обозначения декартовых координат мы использовали буквы в верхнем регистре. Это необходимо, чтобы система могла «различить» случаи, когда они выступают как переменные (приэтом мы будем обозначать их символами в нижнем регистре) и как функции от сферических координат.Функцию ц(р, 6, ф) мы можем представить как сложную, в которой промежуточными аргументами являются р, 9, ф, а независимыми переменными — х, у, z.
Используя приведенные выше формулы, мы найдем частные производные и(р, 9, ф) по декартовым координатам. Повторное дифференцирование полученных функций как сложных даст входящие в выражение лапласианачастные производные второго порядка по х, у и г.Итак, в первую очередь мы должны представить сферические координаты как функции от координат декартовых. Для начала попробуем объединить уравнения, их связывающие, в систему и решить ее с помощью оператора solve относительно р, 9 и ф. Однако полученный при этомрезультат окажется мало приемлемым, так как в нем будет использоваться функция atan2, неимеющая аналогов в «бумажной» математике. Следовательно, выразить р, 6 и ф через х, у и z мыдолжны самостоятельно.Чтобы найти выражение для р, возведем все три уравнения в квадрат и сложим их:X 2 + Y2 + Z simplify -» рТаким образом, г зависит от декартовых координат следующим образом:/222р :=\/х + у + zТак как идентификаторы г как переменной и как функции от х, у, z должны различаться, мы использовали обозначение «р'».Формулу для ф легко найти, разделив выражение для у на выражение для х:Y— simplifyш(ф)sinCOS (ф)Отсюда:ф':=аЬт-Формулу для 9 находим, выразив данную переменную из выражения для Z (учтя полученнуюранее формулу для р):11.2.
Задачи, связанные с вычислением производной * 3 9 1Z9':=acos[—|/ 2 2 2l^x +у +z)Уже известны все величины для того, чтобы можно было найти частные производные и(р, 0, ф)по х, у, г:Г\ Г\Г\ Гах ; Ueду ) \дЬд ,)-Р'Эг )+(д— и\дЪТеперь данные выражения нужно аналитически вычислить. При этом необходимо учесть, чтов выданные в качестве ответа выражения могут входить как сферические координаты, так и декартовы. Последние нужно заменить, используя заданные выше выражения перехода.substitute, х = Х, у = Y , z = Zassume, p > 0,9 = RealRange(O, я) ->simplify+ 2-sin(0) -cos^) + СОБ(Ф) -COS(0)р-П -substitute ,x = X,y=Y,z=ZV=U3assume, p > 0,0 = RealRange(0, я) ->simplify22>-sin^)p-(-l + 2-sin(0) + cos(e) )substitute ,x = X,y=Y,z=Zuz := u zassume,p > O,0 = RealRange(O,rc) ->simplify/А/0ч -1 + 2-sin(9)2 + cos(e) 2sin(0)Чтобы максимально упростить результат аналитического расчета, следует задействовать оператор simplify.
Однако по умолчанию при этом в ответы войдут функции знака выражения — csgn(знак комплексного выражения) и signum (знак действительного выражения). Это происходитиз-за того, что в ряде алгебраических преобразований знак должен учитываться. Чтобы убратьданные функции, символьному процессору нужно показать, в каких пределах изменяются переменные. Сделать это можно с помощью оператора assume. Анализ выражений ответа показывает,что в качестве аргументов функций знаков выступают р и sin(6).
Однако р может быть только392 •Глава 1 1 . Вычисление производныхбольше или равна 0 (см. уравнение, связывающее р и декартовы координаты). То же присущеи sin(G), так как азимутальный угол в случае использованной формы задания сферической системы координат изменяется от 0 до п. Используя операторы неравенства и модификатор RealRange,«сообщаем» эти данные аналитическому процессору. При этом функции знака исчезнут.В той форме, в которой мы получили частные производные функции и(р, 6, ф) по х, у, z, они также являются сложными функциями.
Чтобы вычислить частные производные второго порядка,дифференцируем их, используя уже хорошо знакомые формулы:PЧтобы найти лапласиан, аналитически вычисляем сумму частных производных второго порядка.Так как в выражение ответа могут входить декартовы координаты, их нужно заменить соответствующими им выражениями в сферической системе. Чтобы исключить из ответа функции знака, показываем символьному процессору области изменения г и q.substitute ,x = X,y = Y,z = Zu"x+u"ysimzPlify-»Фсов(ф>мп(в)assume, p > 0,0 = RealRange(o, n)simplifyОбратите внимание, что при расчете последнего выражения оператор simplify был задействовандважды.
Это довольно эффективный ход — проводить упрощение не один раз, а после каждойоперации — поэтому его нужно запомнить.Полученный результат весьма компактен. Однако верен ли он? Чтобы проверить это, воспользуемсяготовым выражением оператора Лапласа в сферической системе координат:AДи :=1 526p —и1+2psin(e)2Ли simplify -> 4-cos (•)>)• sin(8)Ответы сошлись. Следовательно, задача была решена верно.11.3.
Численные методыдифференцированияВ отличие от численного решения уравнений, о котором мы говорили в гл. 8, алгоритм численного дифференцирования, в общем-то, совсем не сложно построить и самостоятельно.Для этого достаточно вспомнить знакомую со школы формулу определения производной:11.3. Численные методы дифференцирования• 393f(x+h)-f(h)— f(x) = limdxh-> 0Очевидно, что, выбрав точку, производная в которой должна быть вычислена, и определившись с величиной приращения, по данной формуле вполне можно вести расчет.Столь же очевидно, что чем меньше величина h, тем точнее будет приближение производной. Впрочем, последнее утверждение не лишним будет и проверить.
Для этого построим таблицу влияния размера шага на точность вычисления производной (табл. 11.1).Для того же, чтобы можно было контролировать величину ошибки, будем дифференцировать функцию в точке, значение производной в которой хорошо известно. Выберем для этого, например, синус в точке п. Точная величина производной для этого случая хорошо известна: —1.f(x) := sin(x)D(a,h):=f(a + h) - f(a)Таблица 1 1 . 1 . Влияние величины приращения на точность численного дифференцированияh1 D(n,h)|-I|0.1F(7i+h)-0.099833416646828D(ii,h)-0.9983341664682820.01-9.99983333416633Е-3-0.9999833334166451.66665833547519E-50.001-9.99999833333109Е-4-0.9999998333332311.66665833547519E-51.66583353171768E-30.0001-9.99999998334219Е-5-0.9999999983354441.66665833547519E-50.00001-9.99999999977639Е-6-0.9999999999898841.01155750442672E-110.000001-1.00000000001715Е-6-1.000000000139611.39611433525033E-100.0000001-9.99999997138813Е-8-0.999999998363421.63658053775606E-90.00000001-9.99999981676465Е-9-0.9999999939225296.07747097092215E-90.000000001-9.99999960279736Е-10-1.000000082740378.27403709990904E-80.0000000000000010-0.8881784197001250.111821580299875Анализ полученной таблицы дает довольно неожиданные результаты: оказывается,точность численного дифференцирования совершенно не обязательно тем выше, чемменьше выбранное приращение.
При шаге 10~15 ошибка более чем в 100 раз выше, чемпри максимальном приращении 10"'. Максимум же точности получается при шаге величиной 10~5. В чем же причина такой странной зависимости точности от приращения?Действительно, в идеальном случае, чем меньше приращение h, тем точнее приближается производная. Однако при использовании компьютера возникает проблема,связанная с тем, что величины и функции вычисляются не абсолютно точно, а существует некоторая расчетная ошибка.
Чем меньшее значение принимает функция, темсильнее сказывается на относительной точности погрешность численных расчетов прималых h разность f(x+h)-f(x) столь невелика, что значение производной будет зависеть от округления весьма заметно. При h, близких по порядкам к точности округления системы, значение производной будет определяться, прежде всего, именно тем, чтоотносительная ошибка при вычислении разности f(x+h)-f(x) сопоставима с самойразностью.
Из-за этого и возникает такая на первый взгляд парадоксальная ситуация3 9 4 •:• Глава 1 1 . Вычисление производныхзначительной неточности определения производной при минимально возможном приращении.В теории численных методов доказывается, что значение приращения, при которомдостигается максимальное приближение к производной, является функцией порогаокругления системы. Так, для используемой нами формулы предела отношений при1/3ращений Ь опт =е , где е — порог точности округления. Так как е в Mathcad равняется17510~ , то, очевидно, дифференцирование следует проводить при h=10~ . Вывод этот полностью подтверждается полученной нами ранее таблицей. При h=l О"5 ошибка начинается только с 11-го знака после запятой, что является очень и очень неплохим приближением.Для иллюстрации вышесказанного попробуем построить график зависимости величины ошибки от значения приращения.
Так как и приращения, и ошибки различаются намного порядков, график следует строить в логарифмическом режиме по обеим шкалам(рис. 11.5).D(a,h) :=i:=0.. 1000f(a + h) - fl»h.:=10-15+1.510 2-iРис. 11.5. График зависимости величины ошибки от значения приращенияДанный график подтверждает наличие максимума точности в окрестности точки приращения 10~5. Помимо того, из него видно, насколько хаотично распределены точкизначений ошибки при малых h. Данная хаотичность объясняется влиянием ошибки округления, имеющей, до определенной степени, случайный характер.Существует довольно значительное количество формул, позволяющих более точнои легко, чем по методу конечных разностей, выполнять численное дифференцирование.