2016-06-14 79 views
-4

我已方程Ax = b的有下列11×11個線性系統:求解線性方程組的病態系統LAPACK&CO

A = { 
{1.0000000000000000, 8.0000000000000000, 6.0000000000000000, 12.0000000000000000, 24.0000000000000000, 24.0000000000000000, 8.0000000000000000, 6.0000000000000000, 24.0000000000000000, 24.0000000000000000, 24.0000000000000000}, 
{4.5999999999999996, 41.8531411531233601, 33.0479488942856037, 87.8349057232554173, 149.3783917109033439, 195.3689938163366833, 121.0451669808013690, 48.8422484540841708, 223.6406089026404516, 851.8470736603384239, 269.3015780207464900}, 
{21.1599999999999966, 218.9606780479085160, 182.0278210198854936, 642.9142219510971472, 929.7459962556697519, 1590.3768227003254196, 1831.4915561762611560, 397.5942056750813549, 2083.9634145976574473, 30235.1432043200838962, 3021.8058301860087340}, 
{97.3359999999999701, 1145.5240206653393216, 1002.6076877338904296, 4705.8591727678940515, 5786.8317341801457587, 12946.2633183243797248, 27711.6501551604087581, 3236.5658295810949312, 19419.1186238102454809, 1073154.9275125553831458, 33907.3782725576675148}, 
{447.7455999999998539, 5992.9723163999815370, 5522.3546042079124163, 34444.8913989153879811, 36017.8173980603314703, 105387.4349242659372976, 419295.1650431178859435, 26346.8587310664843244, 180954.3130575636751018, 38090161.8577392920851707, 380471.2698060897528194}, 
{0.0000000000000000, 34.2801357124991952, 168.4702728821191613, 2101.6181209908259007, 1236.1435394200643714, 6631.0420254749351443, 38374.2674650820554234, 4069.0485156323466072, 28291.8793721561523853, 7044717.1197200166061521, 60211.4334496619121637}, 
{2059.6297599999993508, 31353.0895356311411888, 30417.0821226643129194, 252121.9823892920394428, 224178.4848274685500655, 857893.2134182706940919, 6344206.6583608603104949, 214473.3033545676735230, 1686197.1981563565786928, 1351958038.0734937191009521, 4269229.7229307144880295}, 
{0.0000000000000000, 179.3414198404317403, 927.9328280691040618, 15382.9524602928686363, 7693.8805767663707229, 53979.1670196200575447, 580627.4516345988959074, 33123.5797620395824197, 263633.8804078772664070, 250042569.2999326586723328, 675626.4184535464737564}, 
{0.0000000000000000, 938.2502198978935439, 5111.0461132262771571, 112596.6815912620077142, 47887.4794405465727323, 439410.6478194649680518, 8785268.3545934017747641, 269638.3520710353623144, 2456635.0642409822903574, 8874917956.1941699981689453, 7581135.8600852200761437}, 
{0.0000000000000000, 938.2502198978935439, 0.0000000000000000, 56298.3407956310038571, 23943.7397202732863661, 319571.3802323381532915, 8785268.3545934017747641, 0.0000000000000000, 269630.6777825467870571, 3293783983.7421655654907227, 1735440.7390556528698653}, 
{0.0000000000000000, 70.9608494071368625, 1546.2151390406352220, 34063.2210755480555235, 13279.8613116998949408, 129911.1650312914862297, 2657756.2850107550621033, 183537.2854802548536099, 1654054.3836708476301283, 5487391301.6329326629638672, 5049794.3807012736797333} 
}; 


b = {1, 6.167551546217714, 39.66265463865314, 267.9960092725794, 1918.2310370808632, 137.49061855461255, 14662.396462231256, 1216.4598834815756, 11424.520672986631, 3808.17355766221, 6082.299417407878}; 

矩陣被清楚地病態,雖然正確的解決方案能夠以發現mathematica:

x = {0.0775277, 0.00771443, 0.087553, 0.0208838, 8.47931*1e-7, 0.00197285, 0.0000611365, 0.00187375, 0.000283606, 3.82771*1e-9, 0.000788588}; 

我現在想解決在C程序中使用這個和許多其他類似矩陣的系統。 我已經試過幾乎所有LAPACK函數求解線性方程系統,特別是:

  • dgesv
  • dsgesv
  • dgels
  • dgelss
  • dgelsy

,但他們都給出了嚴重錯誤的結果。

在這一點上,我並不期望從編程的角度來看有任何錯字/錯誤,因爲嘗試使用條件良好的矩陣,我會得到正確的結果。 我想這是概念上的東西,也許我必須使用其他工具..有什麼我可以做的,以找到正確的解決方案與數學圖書館的一些例程?

+2

在這種情況下發布的價值觀是完全沒有必要的,而不是你應該張貼你如何** **嘗試使用* LAPACK *來解決這個問題。另外,如果你知道結果是錯的,它必須是因爲你可以用不同的方法解決這個問題,你爲什麼要在[tag:c]中使用它? –

+1

你在問一個數學問題,而不是一個編程問題。 – kangshiyin

回答

1

解決不良編碼的線性方程組通常很困​​難。至少你不能使用這些一步式LAPACK API來獲得滿意數值錯誤的答案。

作爲一個好的開始,您可以使用截斷SVD方法來獲得更穩定的數值結果。

https://en.m.wikipedia.org/wiki/Linear_least_squares_(mathematics)

這種方法是計算強度最大,但是特別有用的,如果正規方程矩陣,XTX,是非常病態(即,如果其狀態數乘以該機器的相對四捨五入錯誤是相當大的)。在那種情況下,在反演中包括最小的奇異值僅僅爲解決方案增加了數字噪聲。這可以通過截斷SVD方法得到解決,給出一個更加穩定和準確的答案,明確地將所有奇異值設置爲低於某個閾值並忽略它們,這是一個與因子分析密切相關的過程。

更有效的方法可能涉及通過找到一個預調節矩陣求解之前的矩陣充分條件。您需要對原始矩陣的結構有一些瞭解。在下面的討論中你可以找到更多的想法。

https://www.researchgate.net/post/How_can_I_solve_an_ill-conditioned_linear_system_of_equations