Программирование: Несоответствие типов в Maple

[imperfect]

Прошу помощи с такой ошибкой, т.к. хелпов начитался до посинения а ответа нет. Есть простенькая программа:

Код: Выделить всё

> reset; Digits := 15;
#Computations in atomic units

> Enrg := proc (n::integer, l::integer, R::numeric) if n < 0 then error "Main quantum number is less than zero, %1", n elif n-1 < l then error "Illegal orbital number, %1", l elif R <= 0 then error "Illegal radius, %1", R else return evalf(BesselJZeros(l+1/2, n)^2/R^2) end if end proc;
 
> Func := proc (n::integer, l::integer, R::numeric, x) if n < 0 then error "Main quantum number is less than zero, %1", n elif n-1 < l then error "Illegal orbital number, %1", l elif R <= 0 then error "Illegal radius, %1", R else return evalf(sqrt(2)*BesselJ(l+1/2, BesselJZeros(l+1/2, n)*x/R)/(R*BesselJ(l+3/2, BesselJZeros(l+1/2, n))*sqrt(x))) end if end proc;

> StrDipole := proc (n1::integer, l1::integer, n2::integer, l2::integer, R::numeric) local Dipole; if 2 <= abs(l2-l1) then error "It's not a dipole transition, l1 = %1, l2 = %2", l1, l2 else Dipole := evalf(int(Func(n1, l1, R, x)*x^3*Func(n2, l2, R, x), x = 0 .. R))^2; return (2/3)*(Enrg(n2, l2, R)-Enrg(n1, l1, R))*max(l1, l2)*Dipole(n1, l1, n2, l2, R)/(2*l1+1) end if end proc;


Программа состоит из 3х функций, использующих в качестве входные параметры типа integer и numeric. Целью является получение значений последней функции StrDipole в зависимости от заданного последнего параметра. Если написать, например, StrDipole(1,0,2,1,5) нужное число получаем. Но при попытке построить график

Код: Выделить всё

> plot(StrDipole(1, 0, 2, 1, x), x = 1 .. 20);
Error, invalid input: StrDipole expects its 5th argument, R, to be of type numeric, but received x


Как быть?

[xyz]

> plots[listplot]([seq(StrDipole(1, 0, 2, 1, p), p = 1 .. 20)]);
Проработайте вопрос организации процедур в Мапл. Нет времени конкретизировать сказанное. В Мапл нет reset, но есть restart.

[imperfect]

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

[xyz]

Maple имеет достаточно развитый встроенный 4GL язык и хорошо приспособлен для программирования различного типа задач. В справке программирование не описано, с ним лучше всего познакомиться в книгах Аладьева. Я не знаю лучше книг по программированию в Maple, чем книги этого автора и это не только мое мнение. Там вам все станет ясно.

[xyz]

> h:=0.01: plots[listplot]([seq(StrDipole(1, 0, 2, 1, p), p = [(1+k*h)$k=0..4/h])]);

[xyz]

Кстати, определения процедур также следует кодировать верно, например,

Dipole := (n1,n2,l1,l2,R)-> evalf(int(Func(n1, l1, R, x)*x^3*Func(n2, l2, R, x), x = 0 .. R))^2;

Разберитесь также с реализованной в процедурах логикой. Вообще говоря, такиго типа алгоритмы легко и эффективнее реализуются единой процедурой. Дальше уж сами, пожалуйста.

[mois]

Если Вы хотите вывести командой plot график сложной процедуры, то возможны следующие два правильных способа.

1) задать для plot процедуру не как выражение, а как процедуру и указать диапазон изменения x как 1..20, а не как x=1..20:

plot(x->StrDipole(1, 0, 2, 1, x), 1 .. 20);

2) задать для plot процедуру как выражение, но обязательно задержать вычисление параметров процедуры с помощью одинарных кавычек ' (ибо Maple пытается вычислить значения параметров до задания им значений функцией plot) и указать диапазон изменения x как x=1..20:

plot('StrDipole(1, 0, 2, 1, x)', x = 1 .. 20);

или

plot('StrDipole'(1, 0, 2, 1, x), x = 1 .. 20);

[imperfect]

Большое спасибо за ответы, все заработало! Насколько я понял из чтения документации, причина таких ошибок - вычисление maple параметров до присвоения им значений (кавычки этому препятствуют, я сам пытался через модификатор uneval, но не получилось), что отличает его от других языков.

[xyz]

Неверное утверждение о единственно возможных способах вывода графиков процедур. Один я уже привел. Имеются и другие. Наиболее типичные см. в примерах

> P:= proc(x,y) sin(x)+cos(y)+sin(x)*cos(y) end proc;
> plot3d(P(x,y), x=0..Pi, y=0..Pi);
> P1:= proc(x) sin(x)+cos(x) end proc;
> plot(P1(x), x=0..Pi);

Детальнее о процедурах и особенностях их использования см. в вышеуказанной ссылке. Cсылка на сложность процедуры несостоятельна, если процедура закодирована корректно.

[mois]

Неверное утверждение о единственно возможных способах вывода графиков процедур. Один я уже привел. Имеются и другие. Наиболее типичные см. в примерах

> P:= proc(x,y) sin(x)+cos(y)+sin(x)*cos(y) end proc;
> plot3d(P(x,y), x=0..Pi, y=0..Pi);
> P1:= proc(x) sin(x)+cos(x) end proc;
> plot(P1(x), x=0..Pi);

Cсылка на сложность процедуры несостоятельна, если процедура закодирована корректно.

Вы привели примеры простейших процедур, состоящих из одного выражения. Для таких простейших процедур команда plot работает без кавычек.

А вот примеры Ваших же процедур, чуть усложненных, в которые добавлен оператор if:

P:= proc(x,y)
if x<y then
sin(x)+cos(y)+sin(x)*cos(y);
else
sin(x)-cos(y)-sin(x)*cos(y);
end if;
end proc;

P1:= proc(x)
if x<2 then
sin(x)+cos(x);
else
sin(x)-cos(x);
end if;
end proc;

Для таких процедур команды

plot3d(P(x,y), x=0..Pi, y=0..Pi);
plot(P1(x), x=0..Pi);

являются ошибочными и будут выдавать ошибку.

А вот по два варианта правильных команд, которые должны работать для процедур любой сложности:

plot3d('P'(x,y), x=0..Pi, y=0..Pi);
plot3d(P, 0..Pi, 0..Pi);

plot('P1'(x), x=0..Pi);
plot(P1,0..Pi);

[xyz]

Вероятно, вы не писали сложных процедур коли таковыми считаете приведенные. Простые и кодируются очень просто, например:

> P:= (x,y) -> `if`(x<y, sin(x)+cos(y)+sin(x)*cos(y), sin(x)-cos(y)-sin(x)*cos(y));
> P1:= (x) -> `if`(x<2, sin(x)+cos(x), sin(x)-cos(x));
> plot3d(P(x,y), x=0..Pi, y=0..Pi);
> plot(P1(x), x=0..Pi);

> map(type, [P,P1], 'procedure');

[true, true]

> map(type, [P,P1], 'function');

[false, false]

Думаю, что не стоит более обмениваться тривиальными примерами. У меня попросту нет на это времени.