Revisions to How to increase the speed of Collocation method to solve System of differential equations

speed up the code with AffineCovariantNewton method.

Source Link

edited May 28, 2021 at 13:48

68.4k
9
174
489

k = 5.;
tf = 3.; 
c[θ1_?NumericQ, θ2_?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2. Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t]θ1'[t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t]          θ2'[t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

Then, as to the implementation of collocation method, NMinimize is undoubtedly a bad choice. Just use FindRoot with the new-in-12.0 "AffineCovariantNewton" method. I've also turned to Chebyshev–Gauss–Lobatto grid to improve the accuracy of the solution:

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[{arg1, arg2} |-> 
                 With[{t = tpoints, θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
                      Rest /@ {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                               myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}], 
               {#, #} &@{ConstantArray[1., npoints - 1]}];, 
               Method -> "AffineCovariantNewton"]; // AbsoluteTiming
(* {0.864978426819, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

k = 5.;
tf = 3.; 
c[θ1_?NumericQ, θ2_?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2. Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

Then, as to the implementation of collocation method, NMinimize is undoubtedly a bad choice. Just use FindRoot. I've also turned to Chebyshev–Gauss–Lobatto grid to improve the accuracy of the solution:

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[{arg1, arg2} |-> 
                 With[{t = tpoints, θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
                      Rest /@ {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                               myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}], 
               {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

k = 5.;
tf = 3.; 
c[θ1_?NumericQ, θ2_?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2. Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{θ1'[t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
               θ2'[t] == -k Sin[θ1[t]]/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

Then, as to the implementation of collocation method, NMinimize is undoubtedly a bad choice. Just use FindRoot with the new-in-12.0 "AffineCovariantNewton" method. I've also turned to Chebyshev–Gauss–Lobatto grid to improve the accuracy of the solution:

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[{arg1, arg2} |-> 
                 With[{t = tpoints, θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
                      Rest /@ {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                               myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}], 
               {#, #} &@{ConstantArray[1., npoints - 1]}, 
               Method -> "AffineCovariantNewton"]; // AbsoluteTiming
(* {0.426819, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

simplify the code a bit.

Source Link

edited May 28, 2021 at 12:46

xzczd ♦

68.4k
9
174
489

k = 5.;
tf = 3.; 
c[(θ1_)c[θ1_?NumericQ, (θ2_)?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2.*Cos[θ2 Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[
    Function[{arg1, arg2}, |-> 
                 With[{t = tpoints, θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
      Rest /@ With[{t = tpoints},            Rest /@ {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                                   myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]]], 
               {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

|-> is introduced in v12.2, if you're not yet in v12.2, use \[Function] instead.

k = 5.;
tf = 3.; 
c[(θ1_)?NumericQ, (θ2_)?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2.*Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[
    Function[{arg1, arg2}, 
     With[{θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
      Rest /@ With[{t = tpoints}, {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                                   myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]], 
    {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

k = 5.;
tf = 3.; 
c[θ1_?NumericQ, θ2_?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2. Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[{arg1, arg2} |-> 
                 With[{t = tpoints, θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
                      Rest /@ {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                               myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}], 
               {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

|-> is introduced in v12.2, if you're not yet in v12.2, use \[Function] instead.

added 344 characters in body

Source Link

edited May 28, 2021 at 12:05

xzczd ♦

68.4k
9
174
489

k = 5.;
tf = 3.; 
c[(θ1_)?NumericQ, (θ2_)?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2.*Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == ((-k)*Sin[θ1[t]] Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[
    Function[{arg1, arg2}, 
     With[{θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
      Rest /@ With[{t = tpoints}, {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                                   myd . θ2 - ((-k)*Sin[θ1] Sin[θ1])/E^(Sin[t]/10.)}]]], 
    {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - ((-k)*Sin[θ1] Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

k = 5.;
tf = 3.; 
c[(θ1_)?NumericQ, (θ2_)?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2.*Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == ((-k)*Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[
    Function[{arg1, arg2}, 
     With[{θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
      Rest /@ With[{t = tpoints}, {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                                   myd . θ2 - ((-k)*Sin[θ1])/E^(Sin[t]/10.)}]]], 
    {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - ((-k)*Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

k = 5.;
tf = 3.; 
c[(θ1_)?NumericQ, (θ2_)?NumericQ] := 
  NIntegrate[(Sin[θ1 + t]^2.*Cos[θ2 + t])/E^t, {t, 0., tf}, 
   Method -> {Automatic, SymbolicProcessing -> 0}];

sol = NDSolve[{Derivative[1][θ1][t] == θ2[t]/E^(Sin[t]/10.) + c[θ1[t], θ2[t]], 
     Derivative[1][θ2][t] == (-k Sin[θ1[t]])/E^(Sin[t]/10.), 
     θ1[0.] == 1., θ2[0.] == 1.}, {θ1, θ2}, {t, 0., tf}]; // AbsoluteTiming
(* {0.0730725, Null} *)

ref = sol[[1, All, -1]] // ListLinePlot

CGLGrid[{xl_, xr_}, n_] := 1/2 (xl + xr + (xl - xr) Cos[(π Range[0, n - 1])/(n - 1)]);
  
npoints = 11;
tpoints = CGLGrid[{0., tf}, npoints];
myd = dMatrixLagrange[tpoints];

SetAttributes[c, Listable];
θ1init = {1.}; θ2init = {1.};
tst = FindRoot[
    Function[{arg1, arg2}, 
     With[{θ1 = θ1init~Join~arg1, θ2 = θ2init~Join~arg2}, 
      Rest /@ With[{t = tpoints}, {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
                                   myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]], 
    {#, #} &@{ConstantArray[1., npoints - 1]}]; // AbsoluteTiming
(* {0.864978, Null} *)

ListPlot[{tpoints, #}\[Transpose] & /@ Join[{θ1init, θ2init}, tst, 2], 
  PlotRange -> All]~Show~ref

I've made use of the hidden syntax of FindRoot, a more basic approach is:

  help[arg1_List, arg2_] := 
   With[{θ1 = {1.}~Join~arg1, θ2 = {1.}~Join~arg2}, 
    Rest /@ With[{t = tpoints}, 
      {myd . θ1 - (θ2/E^(Sin[t]/10.) + c[θ1, θ2]), 
       myd . θ2 - (-k Sin[θ1])/E^(Sin[t]/10.)}]]

  tst = {arg1, arg2} /. 
   FindRoot[help[arg1, arg2], {#, ConstantArray[1., npoints - 1]} & /@ {arg1, arg2}]

NDSolve`FiniteDifferenceDerivative is a better choice to build myd e.g.:

myd = NDSolve`FiniteDifferenceDerivative[1, tpoints, 
  DifferenceOrder -> 2]["DifferentiationMatrix"]

NDSolve`FiniteDifferenceDerivative[1, tpoints, 
     DifferenceOrder -> "Pseudospectral"]["DifferentiationMatrix"] - 
   dMatrixLagrange[tpoints] // Abs // Max
(* 2.66454*10^-14 *)

added 344 characters in body

Source Link

edited May 28, 2021 at 11:59

xzczd ♦

68.4k
9
174
489

Loading

Source Link

created May 28, 2021 at 10:22

xzczd ♦

68.4k
9
174
489

Loading

Stack Exchange Network

Return to Answer