Решение нелинейной системы с помощью fsolve с несколькими входами в Matlab - PullRequest
Новая
       18

Решение нелинейной системы с помощью fsolve с несколькими входами в Matlab

0 голосов
/ 11 мая 2019

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

Это F, rho и A. Сейчас я установил их все в 1 и смог найти решение, так что, по крайней мере, скрипт работает правильно. 

function f = pressXmanifold(x)
    F = [1,1,1];
    rho = 1;
    A = 1;
    f(1) = x(1)-x(3) - ...
        1/2*rho*(x(4)/F(1))^2*...
        A*(x(4)/x(6)*F(3)/F(1))^(-2)*...
        (1+(F(3)/F(1))^2+3*(F(3)/F(1))^2*((x(4)/x(6))^2-(x(4)/x(6))));
    f(2) = x(2)-x(3) - ...
        1/2*rho*(x(5)/F(2))^2*...
        A*(x(5)/x(6)*F(3)/F(2))^(-2)*...
        (1+(F(3)/F(2))^2+3*(F(3)/F(2))^2*((x(5)/x(6))^2-(x(5)/x(6))));
    f(3) = x(1) - x(3);
    f(4) = x(2) - x(3);
    f(5) = x(4) - x(5);
    f(6) = x(4)+x(5) - x(6);
end

Для F Я бы хотел иметь возможность вводить любой массив 1 x 3 в мою функцию. Для rho я также хотел бы ввести значение. Для A это должно было бы измениться с итерациями, основанными на следующих критериях (заданных как оператор if). 

if F(3)/F(1) <= 0.35 && x(3)/x(2) <= 1
    A = 1;
elseif F(3)/F(1) > 0.35 && x(3)/x(2) <= 0.4
    A = 0.9*(1-x(3)/x(2));
elseif F(3)/F(1) > 0.35 && x(3)/x(2) > 0.4
    A = 0.55;
end

Пока я пробовал работать только с первыми двумя параметрами. f = pressXmanifold(x,F,rho) и удаление первых двух строк в моей функции. Но когда я пытаюсь следовать инструкциям на сайте MathWorks 

options = optimoptions('fsolve','Display','none','PlotFcn',@optimplotfirstorderopt);
fun = @pressXmanifold;
x0 = [1,1,1,1,1,1];
x = fsolve(fun,x0,F,rho,options)

Невозможно выполнить присваивание, поскольку индексирование точек не поддерживается для переменных этого типа. Ошибка в createOptionFeedback (строка 33) параметры. (stopTestOptions {k}) = []; Ошибка в prepareOptionsForSolver (строка 57) optionFeedback = createOptionFeedback (options); Ошибка в fsolve (строка 157) [options, optionFeedback] = prepareOptionsForSolver (options, 'fsolve');

1 Ответ

1 голос
/ 12 мая 2019

Перезаписать pressXmanifold, включая F и rho в качестве входных данных 

function f = pressXmanifold(x1, x2, x3, x4, x5, x6, F, rho)
   if F(3)/F(1) <= 0.35 && x3/x2 <= 1
        A = 1;
    elseif F(3)/F(1) > 0.35 && x3/x2 <= 0.4
        A = 0.9*(1-x3/x2);
    elseif F(3)/F(1) > 0.35 && x3/x2 > 0.4
        A = 0.55;
    end
    f(1) = x1-x3 - ...
        1/2*rho*(x4/F(1))^2*...
        A*(x4/x6*F(3)/F(1))^(-2)*...
        (1+(F(3)/F(1))^2+3*(F(3)/F(1))^2*((x4/x6)^2-(x4/x6)));
    f(2) = x2-x3 - ...
        1/2*rho*(x5/F(2))^2*...
        A*(x5/x6*F(3)/F(2))^(-2)*...
        (1+(F(3)/F(2))^2+3*(F(3)/F(2))^2*((x5/x6)^2-(x5/x6)));
    f(3) = x1 - x3;
    f(4) = x2 - x3;
    f(5) = x4 - x5;
    f(6) = x4+x5 - x6;
end

Затем определить F и rho 

rho = 1;
F =[1,1,1];
fun = @(x)pressXmanifold(x(1), x(2), x(3), x(4), x(5), x(6), F, rho);

Наконец 

options = optimoptions('fsolve','Display','none','PlotFcn',...
    @optimplotfirstorderopt);
x0 = [1,1,1,1,1,1];
x = fsolve(fun,x0,options)
...