Выбор нужных решений и вывод из на график

[smorodina82]

Решаю систему: S:=evalf(solve({ ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2)=((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2), xk=sqrt(R^2-yk^2), (xl-xk)/(2*xk)=zl-zk, (xl-xk)/(2*xk)=(yl-yk)/(2*yk), yl=(xl-xa)*(ys-ya)/(xs-xa)+ya, zl=(xl-xa)*(zs-za)/(xs-xa)+za}, [xk,yk,zk,xl,yl,zl])):
> print(S[4,1],S[4,2],S[4,3]); od; LL:=[seq([S[i,1],S[i,2],S[i,3]],i=1..3)]:
Как вывести на график только четвертые решения xk,yk,zk?

[Kitonum]

Решаю систему: S:=evalf(solve({ ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2)=((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2), xk=sqrt(R^2-yk^2), (xl-xk)/(2*xk)=zl-zk, (xl-xk)/(2*xk)=(yl-yk)/(2*yk), yl=(xl-xa)*(ys-ya)/(xs-xa)+ya, zl=(xl-xa)*(zs-za)/(xs-xa)+za}, [xk,yk,zk,xl,yl,zl])):
> print(S[4,1],S[4,2],S[4,3]); od; LL:=[seq([S[i,1],S[i,2],S[i,3]],i=1..3)]:
Как вывести на график только четвертые решения xk,yk,zk?

Укажите значения всех параметров, т.к. в представленном виде с большим числом параметров система не решается вообще. Непонятно также, какой график имеется в виду, и почему именно четвёртые решения? У Вас есть какая-то дополнительная информация о возможных решениях?

[smorodina82]

restart: with(plots):R:=15: >gelob:=plot3d([R*cos(t),R*sin(t),z],z=0..36,t=-Pi..Pi,axes=boxed):
>xs:=40: ys:=0: zs:=20:xa:=30: za:=0: for ya from 1 to 3 do S:=evalf(solve({ ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2)=((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2), xk=sqrt(R^2-yk^2), (xl-xk)/(2*xk)=zl-zk, (xl-xk)/(2*xk)=(yl-yk)/(2*yk), yl=(xl-xa)*(ys-ya)/(xs-xa)+ya, zl=(xl-xa)*(zs-za)/(xs-xa)+za}, [xk,yk,zk,xl,yl,zl])):
> print(S[4,1],S[4,2],S[4,3]); od;
Меня интересуют все значения [xk,yk,zk], которые каждый раз получаются четвертыми при решении системы. Это координаты точек, которые лежат на цилиндре. Как их вывести на график? Надо сформировать массив? Как это сделать?

[Kitonum]

restart: with(plots):R:=15: >gelob:=plot3d([R*cos(t),R*sin(t),z],z=0..36,t=-Pi..Pi,axes=boxed):
>xs:=40: ys:=0: zs:=20:xa:=30: za:=0: for ya from 1 to 3 do S:=evalf(solve({ ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2)=((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2), xk=sqrt(R^2-yk^2), (xl-xk)/(2*xk)=zl-zk, (xl-xk)/(2*xk)=(yl-yk)/(2*yk), yl=(xl-xa)*(ys-ya)/(xs-xa)+ya, zl=(xl-xa)*(zs-za)/(xs-xa)+za}, [xk,yk,zk,xl,yl,zl])):
> print(S[4,1],S[4,2],S[4,3]); od;
Меня интересуют все значения [xk,yk,zk], которые каждый раз получаются четвертыми при решении системы. Это координаты точек, которые лежат на цилиндре. Как их вывести на график? Надо сформировать массив? Как это сделать?

Вот то, что Вы просили:

restart;
with(plots):
R := 15: xs := 40: ys := 0: zs := 20: xa := 30: za := 0:
L := []:
for ya to 3 do
S := evalf(solve({xk = sqrt(R^2-yk^2), yl = (xl-xa)*(ys-ya)/(xs-xa)+ya, zl = (xl-xa)*(zs-za)/(xs-xa)+za, (xl-xk)/(2*xk) = (yl-yk)/(2*yk), (xl-xk)/(2*xk) = zl-zk, ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2) = ((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2)}, [xk, yk, zk, xl, yl, zl])):
L := [op(L), map(rhs, S[4][1 .. 3])]:
end do:
L;
gelob := plot3d([R*cos(t), R*sin(t), z], z = 0 .. max(L[1, 3], L[2, 3], L[3, 3])+5, t = -Pi .. Pi, style = surface, labels = [x, y, z], numpoints = 10000):
Points := pointplot3d(L, color = [red, green, blue], axes = box, symbol = solidsphere, symbolsize = 17):
display(gelob, Points);
display(gelob, Points, view = [min(seq(L[i, 1], i = 1 .. 3))-1 .. max(seq(L[i, 1], i = 1 .. 3))+1, min(seq(L[i, 2], i = 1 .. 3))-1 .. max(seq(L[i, 2], i = 1 .. 3))+1, min(seq(L[i, 3], i = 1 .. 3))-1 .. max(seq(L[i, 3], i = 1 .. 3))+1]);

[[14.99735648, .281599, 6.7307], [14.98946693, .5620323, 6.74618], [14.97645262, .8401590, 6.77201]]


Боюсь только, что точки могут оказаться не те, которые Вы хотите. Я заметил, что в разных версиях Maple выводит решения в разном порядке (я делал в М12). Поэтому Вам надо чётко сформулировать какие точки Вам нужны без ссылки на номер решения.

[smorodina82]

Спасибо. Буду разбираться. Я тоже заметила, что меня не всегда устраивают 4-е по счету решения. Мне нужны те [xk,yk,zk], которые являются действительными числами и zk>0. Как надо налагать эти условия7

[алексей_алексей]

[quote="smorodina82"][/quote]
Как показалось, Вы решаете систему из 6 уравнений, хотя она сводится к 3. Не знаю, в чём суть задачи, но слышал, что выражения сначала надо упрощать…

[Kitonum]

... Мне нужны те [xk,yk,zk], которые являются действительными числами и zk>0. Как надо налагать эти условия7

Ниже - решение поставленной задачи. В Вашем примере на каждом шаге цикла получаются две точки (эти точки на графике одним цветом).

restart; with(plots):
R := 15: xs := 40: ys := 0: zs := 20: xa := 30: za := 0:
L := []: N := []:
for ya to 3 do evalf(solve({xk = sqrt(R^2-yk^2), yl = (xl-xa)*(ys-ya)/(xs-xa)+ya, zl = (xl-xa)*(zs-za)/(xs-xa)+za, (xl-xk)/(2*xk) = (yl-yk)/(2*yk), (xl-xk)/(2*xk) = zl-zk, ((xk-xs)^2+(yk-ys)^2+(zk-zs)^2)/((xk-xa)^2+(yk-ya)^2+(zk-za)^2) = ((xl-xs)^2+(yl-ys)^2+(zl-zs)^2)/((xl-xa)^2+(yl-ya)^2+(zl-za)^2)}, [xk, yk, zk, xl, yl, zl])):
S := select(s->is(0 < rhs(s[3]) and Im(rhs(s[1])) = 0 and Im(rhs(s[2])) = 0), %):
N := [op(N), nops(S)]:
L := [op(L), seq(map(rhs, S[i][1 .. 3]), i = 1 .. nops(S))]:
n := nops(L):
end do:
L:
gelob := plot3d([R*cos(t), R*sin(t), z], z = -10 .. max(seq(L[i, 3], i = 1 .. n))+20, t = -Pi .. Pi, style = surface, labels = [x, y, z], numpoints = 10000):
Points := pointplot3d(L, color = [`$`(red, N[1]), `$`(green, N[2]), `$`(blue, N[3])], axes = box, symbol = solidsphere, symbolsize = 17):
display(gelob, Points);
display(gelob, Points, view = [min(seq(L[i, 1], i = 1 .. n))-2 .. max(seq(L[i, 1], i = 1 .. n))+2, min(seq(L[i, 2], i = 1 .. n))-2 .. max(seq(L[i, 2], i = 1 .. n))+2, min(seq(L[i, 3], i = 1 .. n))-2 .. max(seq(L[i, 3], i = 1 .. n))+2]);

[smorodina82]

Благодарю. И еще момент: если ya меняется не от 1 до 3, а в других пределах, например, от -10 до 10? А то изменила пределы, и сразу ошибка.

[Kitonum]

Благодарю. И еще момент: если ya меняется не от 1 до 3, а в других пределах, например, от -10 до 10? А то изменила пределы, и сразу ошибка.

Далее давайте уж сами! Разбирайтесь в коде.

[smorodina82]

Спасибо. Буду разбираться)