Skip to content

Commit d7930c8

Browse files
oscardssmithoscardssmith
committed
fix inplace and tableau format
1 parent 325cef9 commit d7930c8

File tree

2 files changed

+56
-96
lines changed

2 files changed

+56
-96
lines changed

‎lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_caches.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -729,7 +729,7 @@ function alg_cache(alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2}, u, rate_proto
729729

730730
# Initialize vectors
731731
dense = [zero(rate_prototype) for _ in 1:2]
732-
ks = [zero(rate_prototype) for _ in 1:5]
732+
ks = [zero(rate_prototype) for _ in 1:6]
733733
du = zero(rate_prototype)
734734
du1 = zero(rate_prototype)
735735
du2 = zero(rate_prototype)

‎lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_tableaus.jl

Lines changed: 55 additions & 95 deletions
Original file line numberDiff line numberDiff line change
@@ -237,69 +237,51 @@ function Rodas4Tableau(T, T2)
237237
#BET2P=0.0317D0
238238
#BET3P=0.0635D0
239239
#BET4P=0.3438D0
240-
# Coefficients for the "a" matrix
241240
A = T[
242-
0 0 0 0 0 0
243-
1.544000000000000 0 0 0 0 0
244-
0.9466785280815826 0.2557011698983284 0 0 0 0
245-
3.314825187068521 2.896124015972201 0.9986419139977817 0 0 0
246-
1.221224509226641 6.019134481288629 12.53708332932087 -0.6878860361058950 0 0
247-
1.221224509226641 6.019134481288629 12.53708332932087 -0.6878860361058950 1 0
241+
0 0 0 0 0 0
242+
1.544 0 0 0 0 0
243+
0.9466785280815826 0.2557011698983284 0 0 0 0
244+
3.314825187068521 2.896124015972201 0.9986419139977817 0 0 0
245+
1.221224509226641 6.019134481288629 12.53708332932087 -0.6878860361058950 0 0
246+
1.221224509226641 6.019134481288629 12.53708332932087 -0.6878860361058950 1 0
248247
]
249-
250-
# Coefficients for the "C" matrix
251248
C = T[
252-
0 0 0 0 0
253-
-5.668800000000000 0 0 0 0
254-
-2.430093356833875 -0.2063599157091915 0 0 0
255-
-0.1073529058151375 -9.594562251023355 -20.47028614809616 0 0
256-
7.496443313967647 -10.24680431464352 -33.99990352819905 11.70890893206160 0
257-
8.083246795921522 -7.981132988064893 -31.52159432874371 16.31930543123136 -6.058818238834054
249+
0 0 0 0 0
250+
-5.6688 0 0 0 0
251+
-2.430093356833875 -0.2063599157091915 0 0 0
252+
-0.1073529058151375 -9.594562251023355 -20.47028614809616 0 0
253+
7.496443313967647 -10.24680431464352 -33.99990352819905 11.70890893206160 0
254+
8.083246795921522 -7.981132988064893 -31.52159432874371 16.31930543123136 -6.058818238834054
258255
]
259-
260-
# Coefficients for the "c" vector
261256
c = T2[0, 0.386, 0.21, 0.63, 1, 1, 1]
262-
263-
# Coefficients for the "d" vector
264-
d = T[0.2500000000000000, -0.1043000000000000, 0.1035000000000000, -0.03620000000000023, 0, 0, 0]
265-
266-
# Coefficients for the "h" (kept as is)
257+
d = T[0.25, -0.1043, 0.1035, -0.0362, 0, 0, 0]
267258
H = T[10.12623508344586 -7.487995877610167 -34.80091861555747 -7.992771707568823 1.025137723295662 0
268259
-0.6762803392801253 6.087714651680015 16.43084320892478 24.76722511418386 -6.594389125716872 0]
269-
270260
RodasTableau(A, C, gamma, c, d, H)
271261
end
272262

273263
function Rodas42Tableau(T, T2)
274264
gamma = convert(T, 1 // 4)
275-
276-
# Create matrices for A and C coefficients
277265
A = T[
278-
0.0 0.0 0.0 0.0 0 0
279-
1.4028884 0.0 0.0 0.0 0 0
280-
0.6581213 -1.320936088384301 0.0 0.0 0 0
281-
7.131197445744498 16.02964143958207 -5.561572550509766 0.0 0 0
282-
22.73885722420363 67.38147284535289 -31.21877493038560 0.7285641833203814 0 0
283-
22.73885722420363 67.38147284535289 -31.21877493038560 0.7285641833203814 1 0
266+
0.0 0 0 0 0 0
267+
1.4028884 0 0 0 0 0
268+
0.6581213 -1.320936088384301 0 0 0 0
269+
7.131197445744498 16.02964143958207 -5.561572550509766 0 0 0
270+
22.73885722420363 67.38147284535289 -31.21877493038560 0.7285641833203814 0 0
271+
22.73885722420363 67.38147284535289 -31.21877493038560 0.7285641833203814 1 0
284272
]
285-
286273
C = T[
287-
0.0 0.0 0.0 0.0 0.0
288-
-5.1043536 0.0 0.0 0.0 0.0
289-
-2.899967805418783 4.040399359702244 0.0 0.0 0.0
290-
-32.64449927841361 -99.35311008728094 49.99119122405989 0.0 0.0
291-
-76.46023087151691 -278.5942120829058 153.9294840910643 10.97101866258358 0.0
292-
-76.29701586804983 -294.2795630511232 162.0029695867566 23.65166903095270 -7.652977706771382
274+
0 0 0 0 0
275+
-5.1043536 0 0 0 0
276+
-2.899967805418783 4.040399359702244 0 0 0
277+
-32.64449927841361 -99.35311008728094 49.99119122405989 0 0
278+
-76.46023087151691 -278.5942120829058 153.9294840910643 10.97101866258358 0
279+
-76.29701586804983 -294.2795630511232 162.0029695867566 23.65166903095270 -7.652977706771382
293280
]
294-
295-
# Create vectors for c and d coefficients
296-
c = T2[0, 0.3507221, 0.2557041, 0.6817790, 1, 1, 1]
297-
d = T[0.2500000000000000, -0.06902209999999998, -0.0009671999999999459, -0.08797900000000025, 0, 0, 0]
298-
299-
# Coefficients for the "H" matrix
281+
c = T2[0, 0.3507221, 0.2557041, 0.681779, 1, 1, 1]
282+
d = T[0.25, -0.0690221, -0.0009672, -0.087979, 0, 0, 0]
300283
H = T[-38.71940424117216 -135.8025833007622 64.51068857505875 -4.192663174613162 -2.531932050335060 0
301284
-14.99268484949843 -76.30242396627033 58.65928432851416 16.61359034616402 -0.6758691794084156 0]
302-
303285
RodasTableau(A, C, gamma, c, d, H)
304286
end
305287

@@ -308,73 +290,51 @@ function Rodas4PTableau(T, T2)
308290
#BET2P=0.D0
309291
#BET3P=c3*c3*(c3/6.d0-GAMMA/2.d0)/(GAMMA*GAMMA)
310292
#BET4P=0.3438D0
311-
312-
# Coefficients for the "A" matrix
313293
A = T[
314-
0.0 0 0.0 0.0 0 0
315-
3.0 0 0.0 0.0 0 0
316-
1.831036793486759 0.4955183967433795 0.0 0.0 0 0
317-
2.304376582692669 -0.05249275245743001 -1.176798761832782 0.0 0 0
318-
-7.170454962423024 -4.741636671481785 -16.31002631330971 -1.062004044111401 0 0
319-
-7.170454962423024 -4.741636671481785 -16.31002631330971 -1.062004044111401 1 0
294+
0 0 0 0 0 0
295+
3 0 0 0 0 0
296+
1.831036793486759 0.4955183967433795 0 0 0 0
297+
2.304376582692669 -0.05249275245743001 -1.176798761832782 0 0 0
298+
-7.170454962423024 -4.741636671481785 -16.31002631330971 -1.062004044111401 0 0
299+
-7.170454962423024 -4.741636671481785 -16.31002631330971 -1.062004044111401 1 0
320300
]
321-
322-
# Coefficients for the "C" matrix
323301
C = T[
324-
0.0 0.0 0.0 0.0 0.0;
325-
-12.0 0.0 0.0 0.0 0.0;
326-
-8.791795173947035 -2.207865586973518 0.0 0.0 0.0;
327-
10.81793056857153 6.780270611428266 19.53485944642410 0.0 0.0;
328-
34.19095006749676 15.49671153725963 54.74760875964130 14.16005392148534 0.0;
329-
34.62605830930532 15.30084976114473 56.99955578662667 18.40807009793095 -5.714285714285717
302+
0 0 0 0 0
303+
-12 0 0 0 0
304+
-8.791795173947035 -2.207865586973518 0 0 0
305+
10.81793056857153 6.780270611428266 19.53485944642410 0 0
306+
34.19095006749676 15.49671153725963 54.74760875964130 14.16005392148534 0
307+
34.62605830930532 15.30084976114473 56.99955578662667 18.40807009793095 -5.714285714285717
330308
]
331-
332-
# Coefficients for the "c" vector
333-
c = T2[0, 3 * gamma, 0.21, 0.63, 1, 1, 1]
334-
335-
# Coefficients for the "d" vector
336-
d = T[0.2500000000000000, -0.5000000000000000, -0.0235040000000000, -0.0362000000000000, 0, 0, 0]
337-
338-
# Coefficients for the "H" matrix
309+
c = T2[0, 0.75, 0.21, 0.63, 1, 1, 1]
310+
d = T[0.25, -0.5, -0.023504, -0.0362, 0, 0, 0]
339311
H = T[25.09876703708589 11.62013104361867 28.49148307714626 -5.664021568594133 0 0
340312
1.638054557396973 -0.7373619806678748 8.477918219238990 15.99253148779520 -1.882352941176471 0]
341-
342313
RodasTableau(A, C, gamma, c, d, H)
343314
end
344315

345316
function Rodas4P2Tableau(T, T2)
346317
gamma = convert(T, 1 // 4)
347-
348-
# Coefficients for the "A" matrix
349318
A = T[
350-
0.0 0.0 0.0 0 0 0
351-
3.0 0.0 0.0 0 0 0
352-
0.906377755268814 -0.189707390391685 0.0 0 0 0
353-
3.758617027739064 1.161741776019525 -0.849258085312803 0 0 0
354-
7.089566927282776 4.573591406461604 -8.423496976860259 -0.959280113459775 0 0
355-
7.089566927282776 4.573591406461604 -8.423496976860259 -0.959280113459775 1 0
319+
0 0 0 0 0 0
320+
3 0 0 0 0 0
321+
0.906377755268814 -0.189707390391685 0 0 0 0
322+
3.758617027739064 1.161741776019525 -0.849258085312803 0 0 0
323+
7.089566927282776 4.573591406461604 -8.423496976860259 -0.959280113459775 0 0
324+
7.089566927282776 4.573591406461604 -8.423496976860259 -0.959280113459775 1 0
356325
]
357-
358-
# Coefficients for the "C" matrix
359326
C = T[
360-
0.0 0.0 0.0 0.0 0.0;
361-
-12.0 0.0 0.0 0.0 0.0;
362-
-6.354581592719008 0.338972550544623 0.0 0.0 0.0;
363-
-8.575016317114033 -7.606483992117508 12.224997650124820 0.0 0.0;
364-
-5.888975457523102 -8.157396617841821 24.805546872612922 12.790401512796979 0.0;
365-
-4.408651676063871 -6.692003137674639 24.625568527593117 16.627521966636085 -5.714285714285718
327+
0 0 0 0 0
328+
-12 0 0 0 0
329+
-6.354581592719008 0.338972550544623 0 0 0
330+
-8.575016317114033 -7.606483992117508 12.224997650124820 0 0
331+
-5.888975457523102 -8.157396617841821 24.805546872612922 12.790401512796979 0
332+
-4.408651676063871 -6.692003137674639 24.625568527593117 16.627521966636085 -5.714285714285718
366333
]
367-
368-
# Coefficients for the "c" vector
369-
c = T2[0, 0.750000000000000, 0.321448134013046, 0.519745732277726, 1, 1, 1]
370-
371-
# Coefficients for the "d" vector
372-
d = T[0.250000000000000, -0.500000000000000, -0.189532918363016, 0.085612108792769, 0, 0, 0]
373-
374-
# Coefficients for the "H" matrix
334+
c = T2[0, 0.75, 0.321448134013046, 0.519745732277726, 1, 1, 1]
335+
d = T[0.25, -0.5, -0.189532918363016, 0.085612108792769, 0, 0, 0]
375336
H = [-5.323528268423303 -10.042123754867493 17.175254928256965 -5.079931171878093 -0.016185991706112 0
376337
6.984505741529879 6.914061169603662 -0.849178943070653 18.104410789349338 -3.516963011559032 0]
377-
378338
RodasTableau(A, C, gamma, c, d, H)
379339
end
380340

0 commit comments

Comments
 (0)