Fortran Coder

查看: 7926|回复: 4
打印 上一主题 下一主题

[数学库] 关于mkl中lapack的线性代数问题

[复制链接]

13

帖子

6

主题

0

精华

熟手

F 币
127 元
贡献
70 点
跳转到指定楼层
楼主
发表于 2019-5-9 03:03:16 | 只看该作者 |只看大图 回帖奖励 |倒序浏览 |阅读模式
本帖最后由 frankgty 于 2019-9-28 00:10 编辑

想解决的是一个Ax=b的问题,第一次使用mkl,选用的函数是SYSV(不确定是否选择正确)。A(Vollmatrix)是一个2400x2400的对称稀疏矩阵,x则是待验证向量Results,右侧Right hand site为RHS。下面是自己试着编的代码,想先试着输出x,然后再跟Results比较。在红色行触发断点,错误0警告0。因为对函数不太了解,输入参数也不知道是否正确,想请教大神帮看一眼,万分感谢。

  Program kontrolle

  use f95_precision
  use lapack95

  Implicit None

  Integer :: n , nrhs , info
  Integer,allocatable :: ipiv(:)
  Real,allocatable :: A(:,:) , B(:,:)
  n = 2400
  nrhs = 1
  Allocate( A(n,n) , B(n,nrhs) , ipiv(n) )

  Open(12,File="Vollmatix.txt")
  Open(13,File="RHS.txt")

  Read(12,*) A(:,:)
  Read(13,*) B(:,:)

  Call SYSV( A(:,:) , B(:,:) , 'L' , ipiv , info )
  Write(*,*) B(:,:)
  Write(*,*) 'IPIV=', ipiv
  Write(*,*) 'INFO=', info

  End Program kontrolle

RHS.txt

32.81 KB, 下载次数: 0

右手侧b

Results.txt

5.81 KB, 下载次数: 0

x

Vollmatrix.zip

169.98 KB, 下载次数: 0

矩阵A

Vm.zip

145.51 KB, 下载次数: 0

矩阵A

分享到:  微信微信
收藏收藏 点赞点赞 点踩点踩

835

帖子

2

主题

0

精华

大宗师

F 币
3926 元
贡献
2334 点
沙发
发表于 2019-5-9 10:04:48 | 只看该作者
[Fortran] 纯文本查看 复制代码
! MKL库函数求解线性方程组
include 'lapack.f90'
program test 
use lapack95
integer,parameter::n=3
integer i
real(4)::a(n,n),b(n)
integer ipiv(n)
a=reshape([1,0,1,2,1,0,3,3,0],[n,n])
b(:)=[14,11,1]
call sgetrf(n,n,a,n,ipiv,i) !LU分解
call sgetrs( 'N', n, n, a, n, ipiv, b,n, i )
print*,b !1,2,3
pause
end program 

13

帖子

6

主题

0

精华

熟手

F 币
127 元
贡献
70 点
板凳
 楼主| 发表于 2019-9-28 00:06:12 | 只看该作者
本帖最后由 frankgty 于 2019-9-28 00:19 编辑


[Fortran] 纯文本查看 复制代码
include 'lapack.f90'
Program kontrolle

use lapack95

Implicit None

Integer,parameter :: n=2400 , nrhs=1
Integer :: i
Integer,allocatable :: ipiv(:)
Real,allocatable :: A(:,:) , B(:,:)

Allocate( A(n,n) , B(n,nrhs) , ipiv(n) )

Open(12,File="Vm.txt")
Open(13,File="RHS.txt")

Read(12,*) A
Read(13,*) B
Close(12)
Close(13)

Call sgetrf(n,n,A,n,ipiv,i)
Call sgetrs('N',n,n,A,n,ipiv,B,n,i)
Print*,B
Pause
End Program kontrolle



您好,十分感谢您的指导,我参照您的例子改写了代码,出现了错误(0x00007FF93CC2B811 (mkl_core.dll)处(位于 Kontrolle.exe 中)引发的异常: 0xC0000005: 读取位置 0x000001D26396D000 时发生访问冲突。)

835

帖子

2

主题

0

精华

大宗师

F 币
3926 元
贡献
2334 点
地板
发表于 2019-9-28 16:32:11 | 只看该作者
本帖最后由 li913 于 2019-9-28 16:57 编辑

给你的例子,你能得到正确结果不?我自己测试了一下,新版本ivf的结果不对。所以,你需要试试新的代码。
[Fortran] 纯文本查看 复制代码
! MKL库函数求解线性方程组
! 依赖 mkl_lapack95.lib
program test 
use lapack95
integer,parameter::n=3
integer i
real(4)::a(n,n), b(n)
integer ipiv(n)
a=reshape([1,0,1,2,1,0,3,3,0],[n,n])
b(:)=[14,11,1]
call getrf(a,ipiv,i)   !LU分解
call getrs(a,ipiv, b)  !a,b均被覆盖
print*,b !1,2,3
pause
end program 

QQ截图20190928123041.png (114.18 KB, 下载次数: 211)

QQ截图20190928123041.png

13

帖子

6

主题

0

精华

熟手

F 币
127 元
贡献
70 点
5#
 楼主| 发表于 2019-9-28 19:09:23 | 只看该作者
li913 发表于 2019-9-28 16:32
给你的例子,你能得到正确结果不?我自己测试了一下,新版本ivf的结果不对。所以,你需要试试新的代码。[mw ...

非常感谢您的耐心指导,已得到正确结果
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

捐赠本站|Archiver|关于我们 About Us|小黑屋|Fcode ( 京ICP备18005632-2号 )

GMT+8, 2024-11-27 00:04

Powered by Tencent X3.4

© 2013-2024 Tencent

快速回复 返回顶部 返回列表