Блог

У меня есть код pyomo для оптимизации траектории полета самолета. У меня есть реалистичный прогноз ветра от ECMWF. Лучшее решение, которое я могу придумать,-это загрузить данные и построить двумерную интерполяцию Шеннона (для этого я использую mathematica, так как я лучше с ней знаком), как на слайдах 15-18 этой презентации, а затем вставить ее в код Pyomo, который, в свою очередь, выглядит так.

 def Wind_lammda_definition1(model, i):
    return m.Wind_lammda[i,1] ==10.655162991835804*Sinc(6.283185307179586*(-350.   m.lammda[i,1]*180/np.pi))*Sinc(6.981317007977318*(-32.   m.phi[i,1]*180/np.pi))   ...
m.Wind_lammda_const1 = Constraint(m.N, rule = Wind_lammda_definition1)
 

бытие m.lammda и m.phi положение самолета.

Мой вопрос заключается в следующем: есть ли способ уменьшить количество терминов, необходимых для этого двумерного Шеннона? Я использую расширение с тысячами терминов, из которых я показываю только первый. Я хочу знать, могу ли я напрямую интерполировать в python и использовать внутреннюю функцию с положением самолета в pyomo или, если вместо этого я могу, в mathematica использовать более эффективную интерполяцию для вставки в python.

Это код mathematica, который я использую, на всякий случай, если он может быть полезен

 SetDirectory[NotebookDirectory[]];
TX1R = Import["Xnodo_1_of_16.mat"];
stepDivisor1X = 20;
stepDivisor2X = 20;
sinc = Sinc[Pi #] amp; ;
X = Flatten[TX1R[[All, All, 1]]];
Y = Flatten[TX1R[[All, All, 2]]];

xmax = X[[Length[X]]];
xmin = X[[1]];
ymax = Y[[Length[Y]]];
ymin = Y[[1]];
dDeltaXwX = (xmax - xmin)/stepDivisor1X ;(* step *)
dDeltaYwX = (ymax - ymin)/stepDivisor2X;(* step *)
aaa = TX1R //. {x_List} :> x;
tamanio = Dimensions[aaa];
countMax = tamanio[[2]] - 2;
savedX = Table[RandomInteger[i, 5], {i, 1, 5}];
savedX = Table[XposX = aaa[[All, 1]]; YposX = aaa[[All, 2]]; 
   windXVal = aaa[[All, i]];
   windXMat = Transpose[{XposX, YposX, windXVal}];
   ifuncEPSX = Interpolation[windXMat //. {x_List} :> x];
   intDataVectorEPSX = 
    Flatten[Table[{t, u, ifuncEPSX[t, u]}, {t, xmin, xmax, 
       dDeltaXwX}, {u, ymin, ymax, dDeltaYwX}]];
   leEPSX = Length@intDataVectorEPSX;
   vectorXintEPSX = 
    Table[intDataVectorEPSX[[i]], {i, 1, leEPSX - 2, 3}];
   vectorYintEPSX = 
    Table[intDataVectorEPSX[[i]], {i, 2, leEPSX - 1, 3}];
   vectorFintEPSX = Table[intDataVectorEPSX[[i]], {i, 3, leEPSX, 3}];
   interpolatedDataEPSX = 
    Transpose[{vectorXintEPSX, vectorYintEPSX, vectorFintEPSX}];
   shannonInterpolationWindX[t_, u_] = 
    Total[#3*sinc[(t - #1)/dDeltaXwX]*sinc[(u - #2)/dDeltaYwX] amp; @@@ 
      interpolatedDataEPSX], {i, 3, countMax   2}];

For[i = 1, i < lengthLoop   1, i  , 
  savedXTextTemp = ExportString[savedX[[i]], "Text"];
  savedXTextTemp2 = 
   StringReplace[
    savedXTextTemp, {"[" -> "(", "]" -> ")", 
     "t" -> StringJoin["m.lammda[i,1]*180/np.pi"], 
     "u" -> StringJoin["m.phi[i,1]*180/np.pi"]}];
  Print[savedXTextTemp2]; savedXReplaced[[i]] = savedXTextTemp2];

ContourPlot[
 shannonInterpolationWindX[t, u], {t, xmin, xmax}, {u, ymin, ymax}, 
 GridLines -> {{xmax, xmin}, None}, PlotLegends -> Automatic]

ListContourPlot[TX1R //. {x_List} :> x, PlotLegends -> Automatic]