Метод Рунге-Кутта - не получается реализовать для системы.

JustAMan · 14.01.2011, 01:28

Bulinator
Спасибо огромное!! В коде - конечно разберёмся! :-)

Сейчас буду изучать его!)

ewert · 14.01.2011, 01:33

JustAMan в сообщении #399647 писал(а):

Bulinator
Спасибо большое за код! Мне бы сперва порядок вычислений "на пальцах" понять..

Перегрузка -- это дело полезное, но сперва доведите до ума свой вариант, там не так уж и много нуждается в исправлении. Глвное -- замените свои идентификаторы на внятные (в соответствии с обозначениями в условии): на t вместо x1 и на x вместо x2. Переменная x2_t действительно нужна, но обозначение никуда не годится -- замените на что-нибудь осмысленное типа dx. Всё это -- в теле основной программы. Затем вернитесь к коду функции myF -- и Вам сразу должно стать понятно, что он никуда не годится (в т.ч. и набор входных аргументов никуда не годится) и как его следует исправить в соответствии с видом правой части исходной системы. Потом нужно будет выловить ещё пару блох: с уже упомянутой неточностью в схеме и с неверным расположением самого последнего оператора в цикле, но это уже потом.

JustAMan · 17.01.2011, 04:50

ewert в сообщении #399653 писал(а):

JustAMan в сообщении #399647 писал(а):

Bulinator
Спасибо большое за код! Мне бы сперва порядок вычислений "на пальцах" понять..

Перегрузка -- это дело полезное, но сперва доведите до ума свой вариант, там не так уж и много нуждается в исправлении. Глвное -- замените свои идентификаторы на внятные (в соответствии с обозначениями в условии): на t вместо x1 и на x вместо x2. Переменная x2_t действительно нужна, но обозначение никуда не годится -- замените на что-нибудь осмысленное типа dx.

Аа.. спасибо! Так действительно смысле понятнее, чем с моими обозначениями.

Я запустил код Bulinator'a - работает очень хорошо. С "ответом" в книге, действительно сходится)) Разница есть где-то в тысячных - ну это и понятно, формулы-то несколько различаются там. Ну это не столь важно, главное, что код работает верно. Но какого же моё удивление было, когда расширив этот случай до системы 4х уравнений - оно почему-то мне выдало совсем не то, что ожидалось :-)

Тестирую я на системе уравнений. Всё на той же - крепление и к ней подвешены два грузика (между креплениями и грузиками - две пружины, соответственно). Уравнения там такие:
$x_1'' = (-k1*x1 + k2*(x_2-x_1) ) / m1;$
$x_2'' = -k2*(x_2-x_1) / m2;$
k1, k2 - жёсткости пружин. m1, m2 - массы грузиков.

Я через замену представляю эти уравнения четырьмя, но с первой производной и получаю систему:
$x_1'=x_3$
$x_2'=x_4$
$x_3' = (-k1*x1 + k2*(x_2-x_1) ) / m1;$
$x_4' = -k2*(x_2-x_1) / m2;$

В x1 и x2 я получаю ответные значения (как я думаю, что там будет смещение грузиков). И я ожидаю там синусоиду на графике. Но получается только возрастающая прямая линия :lol:

В чём тут может быть ошибка? Тут ведь синусоида должна быть.. хотя может она и есть тут где-то, а я её не вижу :-)

Модифицированный код Bulinator'a теперь такой :)

Код:

struct vectorvalue{

   double x1, x2, x3, x4;
   
    vectorvalue(){
        x1=x2=x3=x4=0;
    }
    vectorvalue(double xx1, double xx2, double xx3, double xx4){
        x1=xx1;
        x2=xx2;
        x3=xx3;
        x4=xx4;
    }
    vectorvalue operator +(vectorvalue a){
        a.x1+=this->x1;
        a.x2+=this->x2;
        a.x3+=this->x3;
        a.x4+=this->x4;
        return a;
    }
    vectorvalue operator -(vectorvalue a){
        a.x1=this->x1-a.x1;
        a.x2=this->x2-a.x2;
        a.x3=this->x1-a.x3;
        a.x4=this->x2-a.x4;
        return a;
    }
    void operator =(vectorvalue a){
        this->x1=a.x1;
        this->x2=a.x2;
        this->x3=a.x3;
        this->x4=a.x4;
    }
    vectorvalue operator *(double a){
        vectorvalue res;
        res.x1=a*this->x1;
        res.x2=a*this->x2;
        res.x3=a*this->x3;
        res.x4=a*this->x4;
        return res;
    }
   
};

vectorvalue function( vectorvalue Y, double m1, double m2, double k1, double k2 ){
   return vectorvalue( Y.x3, Y.x4, ( -k1 * Y.x1 + k2 * ( Y.x2 - Y.x1 ) + forceValue ) / m1, ( -k2 * ( Y.x2 - Y.x1 ) ) / m2 );
}

void calc(){
   
   double harmK1 = 0.1, harmK2 = 0.2, harmM1 = 0.8, harmM2 = 2.2, harmStep = 0.01;
   double harmXi = 0.3, harmXi1 = 1.4, harmYi = 4.0, harmYi1 = 6.8;

   vectorvalue func, k1, k2, k3, k4;
   
   func.x1 = harmXi;
   func.x2 = harmYi;
   func.x3 = harmXi1;
   func.x4 = harmYi1;

   printf( "start: x1=%f, x2=%f, x3=%f, x4=%f\n", func.x1, func.x2, func.x3, func.x4 );

   for( unsigned int i = 0; i < 100; i++ ) {
      
      k1 = function( func, harmM1, harmM2, harmK1, harmK2 );
      k2 = function( func + k1 * ( harmStep / 2 ), harmM1, harmM2, harmK1, harmK2 );
      k3 = function( func + k2 * ( harmStep / 2 ), harmM1, harmM2, harmK1, harmK2 );
      k4 = function( func + k3 * harmStep, harmM1, harmM2, harmK1, harmK2 );      
      func=func+(k1+k2*2+k3*2+k4)*(harmStep/6);
      
      printf( "%i: x1=%f, x2=%f, x3=%f, x4=%f\n", i, func.x1, func.x2, func.x3, func.x4 );
      
   }
   
}

Yu_K · 17.01.2011, 05:43

Странная система уравнений у вас - если подвешены грузы - то должна войти сила тяжести. Если грузы лежат горизонтально и связаны пружинами - то http://dxdy.ru/topic28368.html или http://dxdy.ru/topic37022-90.html.

JustAMan · 17.01.2011, 06:37

Yu_K в сообщении #400990 писал(а):

Странная система уравнений у вас - если подвешены грузы - то должна войти сила тяжести. Если грузы лежат горизонтально и связаны пружинами - то http://dxdy.ru/topic28368.html или http://dxdy.ru/topic37022-90.html.

Ну у меня в принципе такие же уравнения (как в первой ссылке): $Изображение$
:-)

По поводу силы тяжести - я добавлял силу в одно из уравнений:
$mx_1''=-k_1*x_1 + k_2*(x_2-x_1) + F(t)$
Но что-то как была прямая, так и осталось она. А в этом вообще ли дело?

Yu_K · 17.01.2011, 07:57

Как правило уравнения в этой задаче (без учета силы тяжести) должны при сложении и вычитании давать два гармонических уравнения для суммы и разности отклонений - у вас этого нет.

Ну и очень интересно силу тяжести учитываете. :-)

.

Для теста возьмите попроще систему например
http://www.wolframalpha.com/input/?i=y% ... 280%29%3D1

Bulinator · 17.01.2011, 20:12

JustAMan
А что есть forceValue в описании function()?

JustAMan · 17.01.2011, 23:18

Bulinator в сообщении #401231 писал(а):

JustAMan
А что есть forceValue в описании function()?

Забыл удалить :) Это некоторое значение силы, которая влияет на систему в определённый промежуток времени.. (в общем-то на время тестов она нулю у меня равняется - покуда думал, что синусоиду получу в начале и дальше буду экспериментировать с силой).

Bulinator · 17.01.2011, 23:30

JustAMan в сообщении #400987 писал(а):

Но получается только возрастающая прямая линия :lol:

Это не возрастающая линия а маленькая часть синуса. Поставьте значения жесткостей пружин равными не Вашим жалким 0.1 и 0.2 а, скажем, 50 и 60. :-)

Получите свой синус :-)

JustAMan · 18.01.2011, 03:45

Bulinator в сообщении #401306 писал(а):

JustAMan в сообщении #400987 писал(а):

Но получается только возрастающая прямая линия :lol:

Это не возрастающая линия а маленькая часть синуса. Поставьте значения жесткостей пружин равными не Вашим жалким 0.1 и 0.2 а, скажем, 50 и 60. :-)

Получите свой синус :-)

Там почему-то такая странная штука, что от изменения k1 и k2 - все x1 и x2 - почему-то не изменяются)) Такое впечатление, что они и в рассчётах вовсе не учавствуют :-)

Почему-то

Не пойму пока, почему...

Bulinator · 18.01.2011, 10:32

Странно, у меня при значении harmK1=50, harmK2=60 выдается:

(результат1)

Код:

start: x1=0.300000, x2=4.000000, x3=1.400000, x4=6.800000
x1=0.326972, x2=4.062934, x3=3.996161, x4=5.785186
x1=0.379885, x2=4.115687, x3=6.580774, x4=4.765490
x1=0.458406, x2=4.158260, x3=9.110708, x4=3.750733
x1=0.561775, x2=4.190751, x3=11.543622, x4=2.750567
x1=0.688818, x2=4.213352, x3=13.838613, x4=1.774346
x1=0.837958, x2=4.226347, x3=15.956854, x4=0.830980
x1=1.007244, x2=4.230109, x3=17.862192, x4=-0.071192
x1=1.194381, x2=4.225087, x3=19.521704, x4=-0.924522
x1=1.396761, x2=4.211804, x3=20.906210, x4=-1.722165
x1=1.611504, x2=4.190849, x3=21.990723, x4=-2.458176
x1=1.835506, x2=4.162863, x3=22.754837, x4=-3.127597
x1=2.065480, x2=4.128532, x3=23.183046, x4=-3.726523
x1=2.298011, x2=4.088577, x3=23.264984, x4=-4.252154
x1=2.529607, x2=4.043739, x3=22.995587, x4=-4.702836
x1=2.756751, x2=3.994772, x3=22.375182, x4=-5.078074
x1=2.975958, x2=3.942427, x3=21.409476, x4=-5.378535
x1=3.183825, x2=3.887445, x3=20.109481, x4=-5.606031
x1=3.377086, x2=3.830541, x3=18.491345, x4=-5.763482
x1=3.552660, x2=3.772396, x3=16.576111, x4=-5.854867
x1=3.707702, x2=3.713647, x3=14.389394, x4=-5.885152
x1=3.839641, x2=3.654877, x3=11.960997, x4=-5.860211
x1=3.946227, x2=3.596605, x3=9.324455, x4=-5.786725
x1=4.025558, x2=3.539280, x3=6.516525, x4=-5.672075
x1=4.076117, x2=3.483274, x3=3.576624, x4=-5.524221
x1=4.096789, x2=3.428878, x3=0.546233, x4=-5.351573
x1=4.086883, x2=3.376297, x3=-2.531741, x4=-5.162855
x1=4.046141, x2=3.325645, x3=-5.613616, x4=-4.966964
x1=3.974743, x2=3.276951, x3=-8.655602, x4=-4.772826
x1=3.873306, x2=3.230153, x3=-11.614472, x4=-4.589257
x1=3.742871, x2=3.185102, x3=-14.448219, x4=-4.424814
x1=3.584893, x2=3.141566, x3=-17.116695, x4=-4.287665
x1=3.401215, x2=3.099233, x3=-19.582230, x4=-4.185449
x1=3.194041, x2=3.057717, x3=-21.810199, x4=-4.125160
x1=2.965907, x2=3.016569, x3=-23.769556, x4=-4.113032
x1=2.719637, x2=2.975278, x3=-25.433303, x4=-4.154432
x1=2.458303, x2=2.933286, x3=-26.778903, x4=-4.253780
x1=2.185180, x2=2.889997, x3=-27.788620, x4=-4.414469
x1=1.903694, x2=2.844785, x3=-28.449793, x4=-4.638813
x1=1.617372, x2=2.797005, x3=-28.755023, x4=-4.927999
x1=1.329788, x2=2.746008, x3=-28.702284, x4=-5.282071
x1=1.044510, x2=2.691150, x3=-28.294952, x4=-5.699921
x1=0.765043, x2=2.631804, x3=-27.541750, x4=-6.179301
x1=0.494783, x2=2.567370, x3=-26.456612, x4=-6.716854
x1=0.236956, x2=2.497288, x3=-25.058466, x4=-7.308158
x1=-0.005421, x2=2.421046, x3=-23.370940, x4=-7.947795
x1=-0.229590, x2=2.338192, x3=-21.422000, x4=-8.629426
x1=-0.433095, x2=2.248341, x3=-19.243514, x4=-9.345886
x1=-0.613812, x2=2.151184, x3=-16.870764, x4=-10.089289
x1=-0.769989, x2=2.046493, x3=-14.341902, x4=-10.851146
x1=-0.900264, x2=1.934129, x3=-11.697361, x4=-11.622496
x1=-1.003691, x2=1.814043, x3=-8.979236, x4=-12.394032
x1=-1.079747, x2=1.686280, x3=-6.230631, x4=-13.156251
x1=-1.128346, x2=1.550981, x3=-3.494999, x4=-13.899589
x1=-1.149834, x2=1.408382, x3=-0.815476, x4=-14.614570
x1=-1.144984, x2=1.258815, x3=1.765786, x4=-15.291948
x1=-1.114982, x2=1.102698, x3=4.208269, x4=-15.922847
x1=-1.061409, x2=0.940541, x3=6.473704, x4=-16.498896
x1=-0.986217, x2=0.772929, x3=8.526668, x4=-17.012355
x1=-0.891692, x2=0.600526, x3=10.335122, x4=-17.456234
x1=-0.780424, x2=0.424057, x3=11.870907, x4=-17.824398
x1=-0.655263, x2=0.244308, x3=13.110176, x4=-18.111661
x1=-0.519274, x2=0.062108, x3=14.033760, x4=-18.313864
x1=-0.375689, x2=-0.121676, x3=14.627459, x4=-18.427942
x1=-0.227856, x2=-0.306150, x3=14.882268, x4=-18.451967
x1=-0.079187, x2=-0.490412, x3=14.794507, x4=-18.385180
x1=0.066896, x2=-0.673552, x3=14.365882, x4=-18.228000
x1=0.207016, x2=-0.854676, x3=13.603455, x4=-17.982025
x1=0.337892, x2=-1.032906, x3=12.519539, x4=-17.650000
x1=0.456392, x2=-1.207401, x3=11.131506, x4=-17.235783
x1=0.559581, x2=-1.377364, x3=9.461520, x4=-16.744283
x1=0.644770, x2=-1.542049, x3=7.536194, x4=-16.181386
x1=0.709555, x2=-1.700776, x3=5.386184, x4=-15.553872
x1=0.751858, x2=-1.852936, x3=3.045715, x4=-14.869305
x1=0.769958, x2=-1.997999, x3=0.552055, x4=-14.135924
x1=0.762520, x2=-2.135521, x3=-2.055060, x4=-13.362521
x1=0.728617, x2=-2.265146, x3=-4.734035, x4=-12.558302
x1=0.667742, x2=-2.386616, x3=-7.442057, x4=-11.732757
x1=0.579823, x2=-2.499763, x3=-10.135755, x4=-10.895508
x1=0.465219, x2=-2.604519, x3=-12.771868, x4=-10.056175
x1=0.324720, x2=-2.700911, x3=-15.307909, x4=-9.224224
x1=0.159532, x2=-2.789059, x3=-17.702820, x4=-8.408829
x1=-0.028736, x2=-2.869172, x3=-19.917606, x4=-7.618734
x1=-0.238098, x2=-2.941545, x3=-21.915936, x4=-6.862125
x1=-0.466221, x2=-3.006551, x3=-23.664700, x4=-6.146505
x1=-0.710461, x2=-3.064634, x3=-25.134524, x4=-5.478587
x1=-0.967897, x2=-3.116302, x3=-26.300221, x4=-4.864196
x1=-1.235381, x2=-3.162113, x3=-27.141181, x4=-4.308178
x1=-1.509584, x2=-3.202673, x3=-27.641693, x4=-3.814341
x1=-1.787044, x2=-3.238617, x3=-27.791183, x4=-3.385392
x1=-2.064219, x2=-3.270602, x3=-27.584386, x4=-3.022908
x1=-2.337543, x2=-3.299298, x3=-27.021427, x4=-2.727314
x1=-2.603477, x2=-3.325370, x3=-26.107823, x4=-2.497881
x1=-2.858565, x2=-3.349470, x3=-24.854402, x4=-2.332746
x1=-3.099484, x2=-3.372229, x3=-23.277140, x4=-2.228947
x1=-3.323097, x2=-3.394241, x3=-21.396919, x4=-2.182469
x1=-3.526496, x2=-3.416053, x3=-19.239210, x4=-2.188319
x1=-3.707054, x2=-3.438162, x3=-16.833689, x4=-2.240605
x1=-3.862455, x2=-3.460998, x3=-14.213781, x4=-2.332631
x1=-3.990736, x2=-3.484923, x3=-11.416153, x4=-2.457013
x1=-4.090316, x2=-3.510220, x3=-8.480156, x4=-2.605792

тогда как при harmK1=0.1, harmK2=0.2 имеем

(результат2)

Код:

start: x1=0.300000, x2=4.000000, x3=1.400000, x4=6.800000
x1=0.314045, x2=4.067983, x3=1.408934, x4=6.796612
x1=0.328179, x2=4.135932, x3=1.417984, x4=6.793175
x1=0.342405, x2=4.203846, x3=1.427152, x4=6.789689
x1=0.356722, x2=4.271726, x3=1.436435, x4=6.786154
x1=0.371134, x2=4.339569, x3=1.445835, x4=6.782571
x1=0.385640, x2=4.407377, x3=1.455350, x4=6.778939
x1=0.400241, x2=4.475148, x3=1.464979, x4=6.775258
x1=0.414940, x2=4.542882, x3=1.474723, x4=6.771530
x1=0.429736, x2=4.610578, x3=1.484582, x4=6.767753
x1=0.444632, x2=4.678237, x3=1.494553, x4=6.763928
x1=0.459627, x2=4.745857, x3=1.504638, x4=6.760056
x1=0.474725, x2=4.813438, x3=1.514835, x4=6.756135
x1=0.489924, x2=4.880979, x3=1.525144, x4=6.752167
x1=0.505228, x2=4.948481, x3=1.535565, x4=6.748151
x1=0.520636, x2=5.015942, x3=1.546097, x4=6.744088
x1=0.536150, x2=5.083363, x3=1.556740, x4=6.739978
x1=0.551771, x2=5.150742, x3=1.567493, x4=6.735821
x1=0.567500, x2=5.218079, x3=1.578355, x4=6.731617
x1=0.583339, x2=5.285374, x3=1.589327, x4=6.727365
x1=0.599287, x2=5.352626, x3=1.600407, x4=6.723067
x1=0.615347, x2=5.419835, x3=1.611595, x4=6.718723
x1=0.631520, x2=5.487000, x3=1.622891, x4=6.714332
x1=0.647805, x2=5.554121, x3=1.634294, x4=6.709895
x1=0.664206, x2=5.621198, x3=1.645803, x4=6.705411
x1=0.680722, x2=5.688230, x3=1.657418, x4=6.700882
x1=0.697355, x2=5.755216, x3=1.669138, x4=6.696307
x1=0.714105, x2=5.822156, x3=1.680964, x4=6.691686
x1=0.730974, x2=5.889049, x3=1.692893, x4=6.687020
x1=0.747963, x2=5.955896, x3=1.704926, x4=6.682308
x1=0.765073, x2=6.022695, x3=1.717063, x4=6.677551
x1=0.782305, x2=6.089447, x3=1.729302, x4=6.672749
x1=0.799659, x2=6.156150, x3=1.741643, x4=6.667901
x1=0.817138, x2=6.222805, x3=1.754085, x4=6.663010
x1=0.834741, x2=6.289410, x3=1.766628, x4=6.658073
x1=0.852471, x2=6.355966, x3=1.779271, x4=6.653092
x1=0.870327, x2=6.422472, x3=1.792014, x4=6.648067
x1=0.888312, x2=6.488927, x3=1.804856, x4=6.642997
x1=0.906425, x2=6.555331, x3=1.817796, x4=6.637884
x1=0.924668, x2=6.621685, x3=1.830834, x4=6.632727
x1=0.943042, x2=6.687986, x3=1.843969, x4=6.627526
x1=0.961547, x2=6.754235, x3=1.857201, x4=6.622281
x1=0.980186, x2=6.820431, x3=1.870529, x4=6.616994
x1=0.998958, x2=6.886575, x3=1.883951, x4=6.611663
x1=1.017865, x2=6.952664, x3=1.897469, x4=6.606289
x1=1.036908, x2=7.018700, x3=1.911081, x4=6.600872
x1=1.056087, x2=7.084682, x3=1.924786, x4=6.595413
x1=1.075404, x2=7.150608, x3=1.938583, x4=6.589911
x1=1.094859, x2=7.216480, x3=1.952473, x4=6.584367
x1=1.114454, x2=7.282296, x3=1.966454, x4=6.578781
x1=1.134189, x2=7.348055, x3=1.980526, x4=6.573153
x1=1.154065, x2=7.413759, x3=1.994688, x4=6.567483
x1=1.174083, x2=7.479405, x3=2.008939, x4=6.561772
x1=1.194244, x2=7.544994, x3=2.023279, x4=6.556019
x1=1.214549, x2=7.610525, x3=2.037707, x4=6.550225
x1=1.234998, x2=7.675998, x3=2.052222, x4=6.544390
x1=1.255593, x2=7.741413, x3=2.066824, x4=6.538514
x1=1.276335, x2=7.806768, x3=2.081512, x4=6.532598
x1=1.297224, x2=7.872065, x3=2.096285, x4=6.526641
x1=1.318261, x2=7.937301, x3=2.111143, x4=6.520643
x1=1.339447, x2=8.002477, x3=2.126084, x4=6.514606
x1=1.360783, x2=8.067593, x3=2.141109, x4=6.508529
x1=1.382269, x2=8.132648, x3=2.156216, x4=6.502412
x1=1.403907, x2=8.197641, x3=2.171405, x4=6.496256
x1=1.425698, x2=8.262573, x3=2.186675, x4=6.490060
x1=1.447641, x2=8.327442, x3=2.202025, x4=6.483825
x1=1.469739, x2=8.392249, x3=2.217455, x4=6.477551
x1=1.491991, x2=8.456993, x3=2.232963, x4=6.471239
x1=1.514398, x2=8.521674, x3=2.248549, x4=6.464888
x1=1.536962, x2=8.586291, x3=2.264213, x4=6.458498
x1=1.559683, x2=8.650844, x3=2.279953, x4=6.452071
x1=1.582561, x2=8.715332, x3=2.295769, x4=6.445605
x1=1.605598, x2=8.779756, x3=2.311661, x4=6.439102
x1=1.628795, x2=8.844114, x3=2.327626, x4=6.432561
x1=1.652151, x2=8.908407, x3=2.343665, x4=6.425983
x1=1.675668, x2=8.972633, x3=2.359777, x4=6.419368
x1=1.699347, x2=9.036794, x3=2.375960, x4=6.412716
x1=1.723188, x2=9.100888, x3=2.392215, x4=6.406028
x1=1.747191, x2=9.164914, x3=2.408541, x4=6.399302
x1=1.771359, x2=9.228874, x3=2.424936, x4=6.392541
x1=1.795690, x2=9.292765, x3=2.441399, x4=6.385743
x1=1.820187, x2=9.356588, x3=2.457931, x4=6.378910
x1=1.844849, x2=9.420343, x3=2.474531, x4=6.372041
x1=1.869678, x2=9.484029, x3=2.491197, x4=6.365136
x1=1.894673, x2=9.547646, x3=2.507928, x4=6.358197
x1=1.919836, x2=9.611193, x3=2.524724, x4=6.351222
x1=1.945168, x2=9.674670, x3=2.541585, x4=6.344212
x1=1.970668, x2=9.738077, x3=2.558509, x4=6.337168
x1=1.996338, x2=9.801413, x3=2.575495, x4=6.330090
x1=2.022178, x2=9.864679, x3=2.592543, x4=6.322977
x1=2.048189, x2=9.927873, x3=2.609652, x4=6.315831
x1=2.074372, x2=9.990995, x3=2.626821, x4=6.308651
x1=2.100726, x2=10.054046, x3=2.644049, x4=6.301437
x1=2.127253, x2=10.117024, x3=2.661335, x4=6.294190
x1=2.153953, x2=10.179929, x3=2.678679, x4=6.286910
x1=2.180827, x2=10.242762, x3=2.696080, x4=6.279598
x1=2.207875, x2=10.305521, x3=2.713537, x4=6.272252
x1=2.235098, x2=10.368207, x3=2.731048, x4=6.264875
x1=2.262496, x2=10.430819, x3=2.748614, x4=6.257465
x1=2.290070, x2=10.493356, x3=2.766233, x4=6.250023
x1=2.317821, x2=10.555819, x3=2.783905, x4=6.242550

-- Вт янв 18, 2011 13:05:52 --

JustAMan
вопрос программисту. Вот если в описании

Код:

vectorvalue function( vectorvalue Y, double m1, double m2, double k1, double k2 )

поставить & перед аргументами:

Код:

vectorvalue function( vectorvalue &Y, double &m1, double &m2, double &k1, double &k2 )

это ускорит программу? Хоть немножечко...

Yu_K · 18.01.2011, 12:26

Лучше сравнить конечно с аналитическим решением - задача линейная с постоянными коэффициентами. Или возьмите какой-нибудь тест простой с известным готовым решением.

Joker_vD · 18.01.2011, 13:11

(Оффтоп)

Передавать double, int и прочие "маленькие" типы по ссылке бессмысленно — ссылка весит столько же, сколько и указатель + вы будете терять время на разыменование. А вот vectorvalue уже можно и по ссылке, т.к. передать все double через стек FPU уже никак не выйдет.

Bulinator · 19.01.2011, 00:31

(Оффтоп)

Joker_vD в сообщении #401427 писал(а):

Передавать double, int и прочие "маленькие" типы по ссылке бессмысленно — ссылка весит столько же, сколько и указатель + вы будете терять время на разыменование. А вот vectorvalue уже можно и по ссылке, т.к. передать все double через стек FPU уже никак не выйдет.

Спасибо!

JustAMan · 20.01.2011, 04:48

Bulinator
Выделил код в отдельный файл из своей проги - блин, он и правда другие результаты выдаёт :lol:

Видимо у меня какие-то ошибки, пошёл разбираться. Ещё раз спасибо за код! :-)

Научный форум dxdy

Метод Рунге-Кутта - не получается реализовать для системы.