FORTRAN 矩阵相乘的函数function怎么写 必须用 function 不能用 subroutine

FORTRAN 矩阵相乘的函数function怎么写 必须用 function 不能用 subroutine,第1张

Fortran中本身就具有数组相乘的函数,叫matmul(matrix_a,matrix_b)

--------------------------------------------------

必须要自己写的话可以这样,在主函数中加入函数接口interface,专门针对function返回值为数组的情况。

--------------------------------------------------

program main

implicit none

!!!!定义function的使用接口

interface

function mymatmul(m,n,a,b)

integer:: m, n

real:: a(m,n), b(n,m)

real:: mymatmul(m,m)

end function mymatmul

end interface

real::a(1,2),b(2,1)

data a /1,2/

data b /1,2/

write(,)mymatmul(1,2,a,b)

end program main

--------------------------------------------------

function mymatmul(m,n,a,b)

integer:: m, n

real:: a(m,n), b(n,m)

real:: mymatmul(m,m)

integer:: i, j

do i=1,m !行循环

do j=1,m !列循环

mymatmul(i,j)=dot_product( a(i,:), b(:,j) ) !点乘

end do

end do

return

end function mymatmul

下面是实现Gauss-Jordan法实矩阵求逆。

#include

<stdlibh>

#include

<mathh>

#include

<stdioh>

int

brinv(double

a[],

int

n)

{

int

is,js,i,j,k,l,u,v;

double

d,p;

is=malloc(nsizeof(int));

js=malloc(nsizeof(int));

for

(k=0;

k<=n-1;

k++)

{

d=00;

for

(i=k;

i<=n-1;

i++)

for

(j=k;

j<=n-1;

j++)

{

l=in+j;

p=fabs(a[l]);

if

(p>d)

{

d=p;

is[k]=i;

js[k]=j;}

}

if

(d+10==10)

{

free(is);

free(js);

printf("errnot

inv

");

return(0);

}

if

(is[k]!=k)

for

(j=0;

j<=n-1;

j++)

{

u=kn+j;

v=is[k]n+j;

p=a[u];

a[u]=a[v];

a[v]=p;

}

if

(js[k]!=k)

for

(i=0;

i<=n-1;

i++)

{

u=in+k;

v=in+js[k];

p=a[u];

a[u]=a[v];

a[v]=p;

}

l=kn+k;

a[l]=10/a[l];

for

(j=0;

j<=n-1;

j++)

if

(j!=k)

{

u=kn+j;

a[u]=a[u]a[l];}

for

(i=0;

i<=n-1;

i++)

if

(i!=k)

for

(j=0;

j<=n-1;

j++)

if

(j!=k)

{

u=in+j;

a[u]=a[u]-a[in+k]a[kn+j];

}

for

(i=0;

i<=n-1;

i++)

if

(i!=k)

{

u=in+k;

a[u]=-a[u]a[l];}

}

for

(k=n-1;

k>=0;

k--)

{

if

(js[k]!=k)

for

(j=0;

j<=n-1;

j++)

{

u=kn+j;

v=js[k]n+j;

p=a[u];

a[u]=a[v];

a[v]=p;

}

if

(is[k]!=k)

for

(i=0;

i<=n-1;

i++)

{

u=in+k;

v=in+is[k];

p=a[u];

a[u]=a[v];

a[v]=p;

}

}

free(is);

free(js);

return(1);

}

void

brmul(double

a[],

double

b[],int

m,int

n,int

k,double

c[])

{

int

i,j,l,u;

for

(i=0;

i<=m-1;

i++)

for

(j=0;

j<=k-1;

j++)

{

u=ik+j;

c[u]=00;

for

(l=0;

l<=n-1;

l++)

c[u]=c[u]+a[in+l]b[lk+j];

}

return;

}

int

main()

{

int

i,j;

static

double

a[4][4]={

{02368,02471,02568,12671},

{11161,01254,01397,01490},

{01582,11675,01768,01871},

{01968,02071,12168,02271}};

static

double

b[4][4],c[4][4];

for

(i=0;

i<=3;

i++)

for

(j=0;

j<=3;

j++)

b[i][j]=a[i][j];

i=brinv(a,4);

if

(i!=0)

{

printf("MAT

A

IS:

");

for

(i=0;

i<=3;

i++)

{

for

(j=0;

j<=3;

j++)

printf("%137e

",b[i][j]);

printf("

");

}

printf("

");

printf("MAT

A-

IS:

");

for

(i=0;

i<=3;

i++)

{

for

(j=0;

j<=3;

j++)

printf("%137e

",a[i][j]);

printf("

");

}

printf("

");

printf("MAT

AA-

IS:

");

brmul(b,a,4,4,4,c);

for

(i=0;

i<=3;

i++)

{

for

(j=0;

j<=3;

j++)

printf("%137e

",c[i][j]);

printf("

");

}

}

}

参考资料:

C常用算法程序集-徐士良

下面是实现Gauss-Jordan法实矩阵求逆。

#include <stdlibh>

#include <mathh>

#include <stdioh>

int brinv(double a[], int n)

{ int is,js,i,j,k,l,u,v;

double d,p;

is=malloc(nsizeof(int));

js=malloc(nsizeof(int));

for (k=0; k<=n-1; k++)

{ d=00;

for (i=k; i<=n-1; i++)

for (j=k; j<=n-1; j++)

{ l=in+j; p=fabs(a[l]);

if (p>d) { d=p; is[k]=i; js[k]=j;}

}

if (d+10==10)

{ free(is); free(js); printf("errnot inv\n");

return(0);

}

if (is[k]!=k)

for (j=0; j<=n-1; j++)

{ u=kn+j; v=is[k]n+j;

p=a[u]; a[u]=a[v]; a[v]=p;

}

if (js[k]!=k)

for (i=0; i<=n-1; i++)

{ u=in+k; v=in+js[k];

p=a[u]; a[u]=a[v]; a[v]=p;

}

l=kn+k;

a[l]=10/a[l];

for (j=0; j<=n-1; j++)

if (j!=k)

{ u=kn+j; a[u]=a[u]a[l];}

for (i=0; i<=n-1; i++)

if (i!=k)

for (j=0; j<=n-1; j++)

if (j!=k)

{ u=in+j;

a[u]=a[u]-a[in+k]a[kn+j];

}

for (i=0; i<=n-1; i++)

if (i!=k)

{ u=in+k; a[u]=-a[u]a[l];}

}

for (k=n-1; k>=0; k--)

{ if (js[k]!=k)

for (j=0; j<=n-1; j++)

{ u=kn+j; v=js[k]n+j;

p=a[u]; a[u]=a[v]; a[v]=p;

}

if (is[k]!=k)

for (i=0; i<=n-1; i++)

{ u=in+k; v=in+is[k];

p=a[u]; a[u]=a[v]; a[v]=p;

}

}

free(is); free(js);

return(1);

}

void brmul(double a[], double b[],int m,int n,int k,double c[])

{ int i,j,l,u;

for (i=0; i<=m-1; i++)

for (j=0; j<=k-1; j++)

{ u=ik+j; c[u]=00;

for (l=0; l<=n-1; l++)

c[u]=c[u]+a[in+l]b[lk+j];

}

return;

}

int main()

{ int i,j;

static double a[4][4]={ {02368,02471,02568,12671},

{11161,01254,01397,01490},

{01582,11675,01768,01871},

{01968,02071,12168,02271}};

static double b[4][4],c[4][4];

for (i=0; i<=3; i++)

for (j=0; j<=3; j++)

b[i][j]=a[i][j];

i=brinv(a,4);

if (i!=0)

{ printf("MAT A IS:\n");

for (i=0; i<=3; i++)

{ for (j=0; j<=3; j++)

printf("%137e ",b[i][j]);

printf("\n");

}

printf("\n");

printf("MAT A- IS:\n");

for (i=0; i<=3; i++)

{ for (j=0; j<=3; j++)

printf("%137e ",a[i][j]);

printf("\n");

}

printf("\n");

printf("MAT AA- IS:\n");

brmul(b,a,4,4,4,c);

for (i=0; i<=3; i++)

{ for (j=0; j<=3; j++)

printf("%137e ",c[i][j]);

printf("\n");

}

}

}

以上就是关于FORTRAN 矩阵相乘的函数function怎么写 必须用 function 不能用 subroutine全部的内容,包括:FORTRAN 矩阵相乘的函数function怎么写 必须用 function 不能用 subroutine、用c语言怎么编写输入一个矩阵求其逆矩阵的程序、c语言矩阵求逆等相关内容解答,如果想了解更多相关内容,可以关注我们,你们的支持是我们更新的动力!

欢迎分享,转载请注明来源:内存溢出

原文地址:https://54852.com/zz/9282591.html

(0)
打赏 微信扫一扫微信扫一扫 支付宝扫一扫支付宝扫一扫
上一篇 2023-04-26
下一篇2023-04-26

发表评论

登录后才能评论

评论列表(0条)

    保存