Здравствуйте. Столкнулся с проблемой при решении системы дифференциальных уравнений которая представлена на рисунке в приложении. Эта система содержит два интервала. Из за того что это система уравнений второго порядка и я не знаю значений для первых производных на границах (y'(-20) и y'(10)) то это задача с граничными условиями и в теории может быть решена с помощью функции bvp4c. Проблема в том как прописать эту систему уравнений вместе с условием непрерывности в районе точки 0 (снизу изображения). Я попытался сделать следующее (функция для передачи в bvp4c):



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

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

Вообще это задача посути есть моделирование электрического потенциала под затвором в МДП транзисторе. Участок где вторая производная равна нулю есть диэлектрик и там действительно согласно закону потенциал меняется линейно. В другой же части моделируется полупроводник. Из за аккумуляции носителей заряда возле границы полупроводник-диэлектрик (точка 0) потенциал там меняется не линейно. Условия y(-0)=y(+0) и a*y'(-0)=c*y'(+0) есть условия непрерывности потенциала и вектора электрической индукции соответственно. Большое спасибо, был бы очень благодарен за какой нибудь маленький примерчик

Похоже что матлаб врятли сможет решить это встроенными функциями. Они все расчитаны на то что производная как и сама функция гладкие а следуя условию a*y'(-0)=c*y'(+0) выходит что первая производная кусочно гладкая. Поэтому врятли это можно решить встроенными функциями, скорее всего придется писать свой алгоритм

Мquote name='X-Shadow' date='Jan 15 2014, 11:05' post='1226262']
Похоже что матлаб врятли сможет решить это встроенными функциями. Они все расчитаны на то что производная как и сама функция гладкие а следуя условию a*y'(-0)=c*y'(+0) выходит что первая производная кусочно гладкая. Поэтому врятли это можно решить встроенными функциями, скорее всего придется писать свой алгоритм
[/quote]

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

И еще о задаче. Ведь речь идет о распределении потенциала по затвором по глубине канала?
Вроде-бы это рассматривалось. Можно, например, посмотреть google.1024.ru или еще где-нибудь.

Кстати, условие y(-0)=y(+0) в этой системе выполняется автоматически при y=0.

В то же время интересно. А это обязательно надо в МАТЛАБе? Это принципиально?

Я задал этот вопрос на хелпе матлаба и там человек предложил решение http://www.mathworks.com.au/matlabcentral/...ential-equation. Опуская матлабовские функции этот мужик предлагает решить подобную задачу следующим образом (как я понял): разбить систему на 2 уравнения и решать их как будто 2 разных дифура с нач условиями. Значения первых производных не известны поэтому для начала приравнять обе к 0 (касательно правой части отрезка это верно т.к там она и равна нулю). Затем решить 2 уравнения отдельно друг от друга и проверить условие в районе точки 0. Если не сходятся то изменить начальное условие (увеличить или уменьшить начальное значение производной) и повторять пока не будут выполняться условия в районе точки 0. Все бы хорошо кроме одного нюанса - я точно знаю что производная в правом участке y'(10)=0 (ну или около нуля т.к в этом месте полупроводник имеет контакт с металлом и на границе т.к это оммический переход потенциалы должны быть равны а следовательно первая производная от функции потенциала тоже равна 0), это условие из физики. Т.е как бы мы это и не меняем. Но вот получается что второе уравнение никак не зависит от изменения параметров в первом. Фактически мы просто подгоняем решение первого уравнения чтобы оно сходилось с решением второго хотя на самом деле решение второго уравнения должно зависеть от граничного условия в точке х=-20 (в этой точке мы прикладываем потенциал U). Или может я чтото не так понял?

Да.
Я примерно так и набросал, правда в МАТКАДе. (Прикрепляю). Нажмите для просмотра прикрепленного файлаНажмите для просмотра прикрепленного файла
Только сначала я взял при х=-20 задал у=5 (в примере это z0), а потом условно-потолочно взял и задал руками y' (это z1) при х=-20.
Потом запустил как решение задачи Коши, то есть на интервале от -20 до 0 с этими начальными условиями, заданными при х=-20.
При этом использовал уравнение для диэлектрика.
В примере использовал функцию rkfixed (в МАТЛАБ должна быть аналогичная "ode45").
В результате получил распределение в диэлектрике и главное - при х=-0 получил значение y (это Rezult1(N1,1)) и y' (это Rezult1(N1,2)).
Здесь вообще Rezult - это массивы решения по точкам (здесь N1+1=100+1 точка решения).

После того как получил y и y' при х=-0, эти значения задал на границе для полупроводника как начальные.
То есть y(+0)=y(-0) и y'(+0)=(c/a)*y'(-0). И дальше продолжил решение для полупроводника от х=0 до х=10.
Подгонкой вручную (подгонял y' (это z1) при х=-20) получил в самой правой точке y=0.

Конечно, такой способ подгонки решения ненормальный. Но раз я неосторожно написал, что попробую - попробовал.
Все-таки это задача краевая и по-хорошему численно ее решать надо-бы наверное методом конечных разностей.
Но это надо еще повозиться. Но непонятно пока - зачем?
Ведь есть пакеты моделирования типа TCAD или хорошие примеры наработаны в flexPDE (прошлая ссылка).

По поводу условий сшивки решения на границе и на контакте. Примерно вроде бы так.
Электрическая индукция не должно претерпевать разрыв на границе двух сред.
А на контакте (в идеальном металле поле равно 0) касательная составляющая напряженности E конечно равна 0, с нормальной сложнее.
А потенциал, который связан через градиент с напряженностью поля может быть и не равен нулю?

Большое спасибо за помощь!

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

Касательно почему не в пакетах TCAD. Я сначала хочу получить аналитический результат а потом сравнить его с симуляцией.

Сейчас пытаюсь реализовать это в MATLAB но пока как то туго. Собственно код:



А сами функции выглядят так:



Но в итоге решение почемуто не сходится. Матлаб просто зависает пытаясь это решить