[PerlChina] Re: 请问哪位用perl写过解方程组的脚本？

没有人回答我的问题，只好自己参照网上python的程序写了个gauss-jordan消元的程序，将增广矩阵作消元后，得到简化行梯矩阵，参数是增
广矩阵的reference。
好像这个组里能真正帮上忙的不多...

sub gauss_jordan
{
print "gauss_jordan elimination!\n";
my $m = &matrix_row($_[0]); print "there are $m rows in the
matrix!\n";
my $n = &matrix_column($_[0]); print "there are $n columns in
the matrix!\n";
my $k = $m;
$k=$n, if($k>$n);

for(my $i=0;$i<$k;$i++) #转换成上三角矩阵
{
my $maxrow = $i;
for(my $j=($i + 1) ;$j<$m;$j++)
{
if(abs($_[0]->[$j][$i]) > abs($_[0]->[$maxrow][$i]))
{
$maxrow = $j;
}
}
my @temp;
@temp = @{$_[0]->[$i]};
@{$_[0]->[$i]} = @{$_[0]->[$maxrow]};
@{$_[0]->[$maxrow]} = @temp;
#换行成功！接下来将该列的余下的元素全部消为0
if(abs($_[0]->[$i][$i]) > 0 )
{
for(my $j=($i + 1) ;$j<$m;$j++)
{
my $c = $_[0]->[$j][$i]/$_[0]->[$i][$i];
for(my $x =$i;$x<$n;$x++)
{
$_[0]->[$j][$x] -= $_[0]->[$i][$x]*$c;
}
}
}
else
{
print "find a o column!\n";
}
} ##成功转换成上三角矩阵
&output_matrix($_[0]);
print "\n";
##下面将上三角矩阵尽量简化，反向将非对角线元素消为0
for(my $i=$k-1;$i>=0;$i--)
{
if(abs($_[0]->[$i][$i])>0)
{
my $c = $_[0]->[$i][$i];
for(my $j=0; $j<$i;$j++)
{
for(my $x=$n-1;$x>=$i;$x--)
{
$_[0]->[$j][$x] -= $_[0]->[$i][$x]*$_[0]-
>[$j][$i]/$c;
}
}
##以上已经将上三角非对角线元素消为0
##以下将对角线元素消为1
for(my $x=$i;$x<$n;$x++)
{
$_[0]->[$i][$x] /=$c;
}
}
else
{
print "the diagonal element is 0!\n";
}
}
&output_matrix($_[0]);
return($_[0]);
}

# &matrix_row($m)
sub matrix_row
{ my $array = $_[0];
my $k=0;
foreach (@{$array})
{
$k++;
# print "$_!\n";
}
return($k);
}

# &matrix_column($m)
sub matrix_column
{
my $array = $_[0];
my $k=0;
foreach (@{$array->[0]})
{
$k++;
# print "$_!\n";
}
return($k);
}

sub output_matrix
{
my $m = &matrix_row($_[0]); # print "there are $m rows in the
matrix!\n";
my $n = &matrix_column($_[0]); # print "there are $n columns in
the matrix!\n";
for(my $i=0;$i<$m;$i++) #转换成上三角矩阵
{
my $line = join("\t",@{$_[0]->[$i]});
print "$line\n";
}
}







On 7月5日, 下午4时20分, Hydraphenix <hydrapheni...@gmail.com> wrote:
> 要做一些数据的计算，使用math-assistant中的Solve_Det函数，要求系数矩阵必须是正方形矩阵。我的数据中不能保证是正方形矩阵。
> 但计算前可以判断矩阵的秩与未知数的个数相等。在网上搜索了半天，都没有好的思路，希望有朋友能指点下。

--
您收到此邮件是因为您订阅了 Google 网上论坛的"PerlChina Mongers 讨论组"论坛。
要向此网上论坛发帖，请发送电子邮件至 perlchina@googlegroups.com。
要取消订阅此网上论坛，请发送电子邮件至 perlchina+unsubscribe@googlegroups.com。
若有更多问题，请通过 http://groups.google.com/group/perlchina?hl=zh-CN 访问此网上论坛。