Matlab, как мне написать заявление, которое даст мне время на xaxis от y = 0.3 - PullRequest
3 голосов
/ 07 апреля 2011
x=[0:.01:10];
y=(x.^2).*(exp(-x));

plot(x,y), grid
y1=max(y);

xlabel('TIME');
ylabel('CONCENTRATION IN BLOOD');
title('CONCENTRATIN OF SUSTANCE IN BLOOD vs TIME');

fprintf('(a) The maximum concentraion is %f \n',y1)

Это моя программа, и мне сложно написать заявление, которое даст мне время, когда y=0.3, пожалуйста, помогите
, спасибо

Ответы [ 5 ]

4 голосов
/ 08 апреля 2011

Одна из ключевых проблем заключается в том, что на вашем графике есть несколько точек, где y = 0.3. Если вы хотите найти все из них простым способом, вы можете выполнить следующие действия:

  • Вычтите 0,3 из вашего вектора y, чтобы точки, которые вы хотите найти, стали пересечениями с нулями.
  • Найдите индексы в вышеприведенном векторе, где происходит смена знака.
  • Для значений y по обе стороны от пересечения нуля вычислите процент от разницы между ними, при котором лежит значение 0,3. По существу, выполняется линейная интерполяция между двумя точками по обе стороны от пересечения нуля.
  • Используйте вышеуказанный процент, чтобы найти соответствующее значение x для пересечения нуля.

И вот код вместе с графиком, чтобы показать найденные точки:

>> yDesired = 0.3;
>> index = find(diff(sign(y-yDesired)));
>> pctOfDiff = (yDesired-y(index))./(y(index+1)-y(index));
>> xDesired = x(index)+pctOfDiff.*(x(index+1)-x(index))

xDesired =

    0.8291    3.9528

>> plot(x,y);
>> hold on;
>> plot(xDesired,yDesired,'r*')
>> xlabel('x');
>> ylabel('y');

enter image description here

1 голос
/ 07 апреля 2011

Если у вас есть символическая панель инструментов в MATLAB, вы можете сделать следующее

syms x
x=solve('x^2*exp(-x)=y')

x=
  (-2)*lambertw(k, -((-1)^l*y^(1/2))/2)

Здесь lambertw - решение для y=x*exp(x), которое доступно как функция в MATLAB. Теперь вы можете определить функцию как

t=@(y,k,l)(-2)*lambertw(k, -((-1)^l*y^(1/2))/2)

lambertw - это многозначная функция с несколькими ветвями. Переменная k позволяет вам выбрать ответвление решения. Вам нужна главная ветвь, следовательно k=0. l (строчная буква L) просто для выбора подходящего квадратного корня из y. Нам нужен положительный квадратный корень, следовательно l=0. Следовательно, вы можете получить значение t или время для любого значения y, используя функцию.

Итак, используя ваш пример, t(0.3,0,0) дает 0.8291.

EDIT

Я забыл, что есть две ветви решения, которые дают вам реальные результаты (ответ gnovice напомнил мне об этом). Итак, для обоих решений используйте

t(0.3,[0,-1],0)

, что дает 0.8921 и 3.9528.

1 голос
/ 07 апреля 2011

Простой ответ:

find(min(abs(y- 0.3))== abs(y- 0.3))

, дающий

ans = 84

, таким образом

x(84)
ans = 0.83000

Теперь, если вы увеличите разрешение в x, выВы также сможете найти решение ближе к аналитическому.

> x=[0.5:.000001:1]; y=(x.^2).*(exp(-x));
> x(find(min(abs(y- 0.3))== abs(y- 0.3)))
ans =  0.82907

Редактировать :
И конечно, чтобы найти все нули:

> x=[0:.01:10]; y=(x.^2).*(exp(-x));
2> find(abs(y- 0.3)< 1e-3)
ans =
    84   396
> x(find(abs(y- 0.3)< 1e-3))
ans =
   0.83000   3.95000
0 голосов
/ 29 июня 2011

В дополнение к уже опубликованным решениям я добавляю другие способы решения проблемы:

f = @(x) (x.^2).*(exp(-x));             %# function handle
y0 = 0.3;

format long

%# find root of function near s0
x1 = fzero(@(x)f(x)-y0, 1)              %# find solution near x=1
x2 = fzero(@(x)f(x)-y0, 3)              %# find solution near x=3

%# find minimum of function in range [s1,s2]
x1 = fminbnd(@(x)abs(f(x)-y0), 0, 2)    %# find solution in the range x∈[0,2]
x2 = fminbnd(@(x)abs(f(x)-y0), 2, 4)    %# find solution in the range x∈[2,4]
0 голосов
/ 08 апреля 2011

Более простой способ найти индекс (и, следовательно, значение x):

[minDiff, index] = min(abs(y-0.3))

minDiff =

  3.9435e-004


index =

    84

 x(index)

ans =

    0.8300
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...