So I need to find Specific Volumes using the Newton Method from the Van der Waal equation over a change of constant Temperature but variant Pressure.
I have been able to make a list of however many iterations of the altered Van der Waal equation for the root finding method from Pressure Min to Pressure Max (3.9atm to 59.2atm). Further, I am able to select one element/expression from that list and find the root/Specific Volume for that iteration.
However, I cannot seem to figure out how to make a list of all the Specific Volumes, ie solve for the root in more than one expression at a time.
This is my code solving for the 10th element in the list (about 12.2atm):
a = 1.36;
b = .003183;
R = .0820578;
T = 333;
inc = (59.2 - 3.9)/60;
c = Table[ vi - ((((atm) + (a/vi^2)) (vi - b) - (R*T))/((atm) - (a/ vi^2) + (2 a*b/vi^3))), {atm, 3.9, 59.2, inc}]
FindRoot[c[[10]] == 0, {vi, 4}, WorkingPrecision -> 20]
The output:
{-((-27.3252 + (3.9 + 1.36/vi^2) (-0.003183 + vi))/(
3.9 + 0.00865776/vi^3 - 1.36/vi^2)) +
vi,
-((-27.3252 + (4.82167 + 1.36/vi^2) (-0.003183 + vi))/(
4.82167 + 0.00865776/vi^3 - 1.36/vi^2)) +
vi,
-((-27.3252 + (5.74333 + 1.36/vi^2) (-0.003183 + vi))/(
5.74333 + 0.00865776/vi^3 - 1.36/vi^2)) + vi,
AND SO ON.............
Also:
{vi -> 2.1930639624928306932}
Which is a correct answer.
I've tried this fix to solve for more than one element in the list(elements one to ten):
FindRoot[c[[1;;10]] == 0, {vi, 4}, WorkingPrecision -> 20]
but this yields an error:
FindRoot::nveq: The number of equations does not match the number of variables in FindRoot[c[[1;;10]]==0,{vi,4},WorkingPrecision->20]. >>
Any ideas? Thanks!
EDIT: I'm sorry, it's working again for the Van der Waals but neither of the other 2. No matter how I change the starting point the output now is the same, odd number repeated and I don't know how to manipulate either FindMinimum nor Foldlist to change that.
My code is
R = 8.314;
A = 507.28;
B = 0.10476;
a = 0.07132;
b = 0.07235;
c = 660000;
T = 333;
iter = 60;
FX = -P - (A (1 - a/v))/v^2 + (R T (1 - c/(T^3 v)) (B (1 - b/v) + v))/
v^2;
FXX = -((a A)/v^4) + (2 A (1 - a/v))/v^3 + (
R T (1 + (b B)/v^2) (1 - c/(T^3 v)))/v^2 + (
c R (B (1 - b/v) + v))/(T^2 v^4) - (
2 R T (1 - c/(T^3 v)) (B (1 - b/v) + v))/v^3;
FXXX = v - FX/FXX;
nn = Table[FXXX, {P, 400, 6000, iter}];
Ps = Table[P, {P, 400, 6000, iter}]
vs = Rest@FoldList[v /. Last@FindMinimum[#2, {v, #1}] &, 12, nn]
ListPlot[Thread[{Ps, vs}], Joined -> True, PlotRange -> {0, 7}]
The output is
{1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77,
1.3665*10^77, 1.3665*10^77, 1.3665*10^77, 1.3665*10^77}
I'm sorry I basically know nothing haha, but I need this done and I just don't know how to do it. Thanks so much for your help.