谁能帮忙把这段Fortran代码转换成C代码,分不够可以再加
如题,程序如下:
(1)读入数据程序
PROGRAM MAIN
COMMON TYP(10),NO(10),NFO(10),NTO(10),NCO(10),VAL(10)
COMMON/B1/T (29, 30), NOD, NBR, M, N
OPEN (3, FILE=’DD.DAT’, STATUS=’OLD’)
READ (3,*)
READ (3,*) NBR
WRITE (*,*)
WRITE (*,*)
DO 10 I=1, NBR
READ (3,2) TYP(I),NO(I),NFO(I),NTO(I),NCO(I),VAL(I)
WRITE(*,4) TYP(I),NO(I),NFO(I),NTO(I),NCO(I),VAL(I)
NOD=MAX0 (NOD, NFO (I), NTO (I))
10 CONTINUE
2 FORMAT (A2, 4I3, G10.3)
4 FORMAT (5X, A2, 3X, 3(I3, 2X), I3.0, G12.4)
M=NOD+2*NBR
N=M+1
CALL FORMT
CALL GAUSS
CALL OUTPUT
END
(2) 形成表格方程
SUBROUTINE FORMT
COMMON TYP(10),NO(10),NFO(10),NTO(10),NCO(10),VAL(10)
COMMON/B1/T (29, 30), NOD, NBR, M, N
DO 10 J=1, NBR
NF=NFO (J)
NT=NTO (J)
MI=NOD+NBR+J
MJ=NOD+J
NI=MI
NJ=MI
C KCL: Ai=0 and KVL: u=AT*v
IF (NF.NE.0) THEN
T (NF, MI) =1.
T (MJ, NF) =-1.
END IF
IF (NT.NE.0) THEN
T (NT, MI) =-1
T (MJ, NT) =-1.
END IF
T (MJ, MJ) =1.0
C VCR: Mv + Ni = u
IF (TYP (J).EQ.’R’.OR.TYP (J).EQ.’V’) THEN
T (MI, MJ) =1
IF(TYP(J).EQ.’R’) T(NI,NJ)=-VAL(J)
IF(TYP(J).EQ.’V’) T(NI,N)=VAL(J)
END IF
IF (TYP (J).EQ.’G’.OR.TYP (J).EQ.’I’) THEN
IF(TYP(J).EQ.’G’) T(MI,MJ)=-VAL(J)
IF(TYP(J).EQ.’I’) T(MI,N)=VAL(J)
T(NI,NJ)=1
END IF
IF (TYP (J).EQ.’VV’.OR.TYP (J).EQ.’VC’) THEN
IF(TYP(J).EQ.’VV’) T(MI,MJ)=1
IF(TYP(J).EQ.’VC’) T(NI,NJ)=1
T(MI,NOD+NCO(J))=-VAL(J)
END IF
IF (TYP (J).EQ.’CC’.OR.TYP (J).EQ.’CV’) THEN
IF(TYP(J).EQ.’CC’) T(NI,NJ)=1
IF(TYP(J).EQ.’CV’) T(NI,NJ)=1
T(NI,NOD+NBR+NCO(J))=-VAL(J)
END IF
10 CONTINUE
END
(3)调用GAUSS求解方程
SUBROUTINE GAUSS
COMMON / B1 / T(29,30),NOD,NBR,M,N
DO 10 K = 1,M
20 IF(ABS(T(I,K)).GT.ABS(T(L,K))) L=I
IF(ABS(T(L,K)).LT.1.E-30) STOP ’电路无唯一解’
IF(L.NE.K) THEN
DO 30 J = K,N
T1 = T(K,J)
T(K,J) = T(L,J)
30 T(L,J) = T1
END IF
DO 40 I = K,M
C = T(I,K)
IF(C.NE.0.0) THEN
DO 50 J = K,N
T(I,J) = T(I,J)/C
50 IF(I.GT.K) T(I,J) = T(I,J)-T(K,J)
END IF
40 CONTINUE
10 CONTINUE
DO 60 I = M-1, 1,-1
DO 60 J = M, I+1,-1
60 T(I,N) = T(I , N)-T(I ,J)*T(J, N)
END
输出子程序
SUBROUTINE OUTPUT
COMMON TYP (10), NO (10), NFO (10), NTO (10), NCO (10), VAL (10)
COMMMON / B1 / T (29, 30), NOD, NBR, M, N
WRITE (*, 2)
2 FORMAT(/5X,’节点电压’)
WRITE (*, 4) (J, T (J, N)), J = 1, NOD)
4 FORMAT(5X, 4(:,’V’,I2,’=’,G12.4,1X))
WRITE (*, 6)
6 FORMAT( /5X,’支路电压’,9X,’支路电流’,9X’支路功率’)
DO 10 J = 1, NBR
P = T (NOD+J, N)*T (NOD+NBR+J, N)
10 WRITE (*, 12) J, T (NOD+J, N), J T (NOD+NBR+J, N), J, P
12 FORMAT (5X, ‘U’, I2,’=’, G12.4,’I’, I2,’=’, G12.4,’P’, I2,’=’, G12.4)
END
[解决办法]
回楼主
http://www.netlib.org/f2c/f2c.h
http://www.netlib.org/f2c/mswin/f2c.exe.gz
http://www.netlib.org/f2c/libf2c.zip
http://www.netlib.org/f2c/libf2c.zip
是源码,需要自己编译,
NMAKE -f Makefile.vc
[解决办法]
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
/* ------------------------
funtion declarations
----------------------------*/
void format();
void gauss();
void output();
void skip_line(FILE* fp);
int read_line_char(FILE* fp, char* firstToken);
int read_line_int(FILE* fp, int* firstToken);
int read_float(double* value);
int read_int(int* value);
int same_string(const char* str1, const char* str2);
int max(int a, int b);
/*
COMMON TYP(10),NO(10),NFO(10),NTO(10),NCO(10),VAL(10)
COMMON/B1/T (29, 30), NOD, NBR, M, N
*/
char TYP[11][3]={ '\0 '};
int N0[11];
int NF0[11];
int NT0[11];
int NC0[11];
double VAL[11];
double T[30][31];
int NOD;
int NBR;
int M;
int N;
/*----------------------------
reading stuffs
------------------------------*/
static char gBuff[255];
static char gDel[] = ", \ "\ '\n!%\t ";
void skip_line(FILE* fp)
{
int nChar=254;
fgets(gBuff, nChar, fp);
}
/*
read a line, return the first token of the line
*/
int read_line_char(FILE* fp, char* firstToken)
{
int nChar = 254;
char* token = NULL;
if(fgets(gBuff, nChar, fp)) {
token = strtok(gBuff, gDel);
if(token) _snprintf(firstToken, 30, token);
else sprintf(firstToken, " ");
return 1;
}
else {
if(feof(fp)) return EOF;
return 0;
}
}
int read_line_int(FILE* fp, int* firstToken)
{
*firstToken = 0;
int nChar = 254;
char* token = NULL;
if(fgets(gBuff, nChar, fp)) {
token = strtok(gBuff, gDel);
if(token) {
*firstToken = atoi(token) ;
return 1;
} else return 0;
}
else {
if(feof(fp)) return EOF;
return 0;
}
}
int read_float(double* value)
{
char* token = strtok(NULL,gDel);
if(!token) return 0;
*value = atof(token);
return 1;
}
int read_int(int* value)
{
char* token = strtok(NULL,gDel);
if(!token) return 0;
*value = atoi(token);
return 1;
}
/*-------------------------------------------------*/
/* ------------------------------------------------
utilities
---------------------------------------------------*/
int same_string(const char* str1, const char* str2)
{
if (stricmp(str1, str2)==0) return 1;
else return 0;
}
int max(int a, int b)
{
return a> b?a:b;
}
//(1)读入数据程序
int main(int argc, char* argv[])
{
char* buff = (char*) malloc(255*sizeof(char));
FILE* fp = NULL;
fp = fopen( "dd.dat ", "rt ");
if(!fp) {
printf( "din 't find file dd.data! ");
exit(1);
}
skip_line(fp);
read_line_int(fp, &NBR);
printf( "\n\n ");
for(int I=1; I <=NBR; ++I) {
read_line_char(fp, buff);
strcpy(TYP[I], buff);
read_int( &N0[I]);
read_int( &NF0[I]);
read_int( &NT0[I]);
read_int( &NC0[I]);
read_float( &VAL[I]);
printf( "\t %s %d %d %d %d %g\n ", TYP[I],N0[I],NF0[I],NT0[I],NC0[I],VAL[I]);
NOD = max(NOD, max(NF0[I], NT0[I]) );
}
fclose(fp);
free(buff);
M=NOD+2*NBR;
N=M+1;
format();
gauss();
output();
}
//(2) 形成表格方程
void format()
{
int NF, NT, MI, MJ, NI, NJ;
for (int J=1; J <=NBR; ++J) {
NF=NF0[J];
NT=NT0[J];
MI=NOD+NBR+J;
MJ=NOD+J;
NI=MI;
NJ=MI;
// KCL: Ai=0 and KVL: u=AT*v
if (NF!=0) {
T[NF][MI] =1;
T[MJ][NF] =-1;
}
if (NT!=0) {
T[NT][MI] =-1;
T[MJ][NT] =-1;
}
T[MJ][MJ] =1.0;
// VCR: Mv + Ni = u
if ( same_string(TYP[J], "R ") || same_string(TYP[J], "V ") ) {
T[MI][MJ] = 1;
if( same_string(TYP[J], "R ") ) T[NI][NJ]=-VAL[J];
if( same_string(TYP[J], "V ") ) T[NI][N]=VAL[J];
}
if ( same_string(TYP[J], "G ") || same_string(TYP[J], "I ") ) {
if( same_string(TYP[J], "G ") ) T[MI][MJ]=-VAL[J];
if( same_string(TYP[J], "I ") ) T[MI][N]=VAL[J];
T[NI][NJ]=1;
}
if (same_string(TYP[J], "VV ")|| same_string(TYP[J], "VC ") ) {
if(same_string(TYP[J], "VV " ) ) T[MI][MJ]=1;
if(same_string(TYP[J], "VC ") ) T[NI][NJ]=1;
T[MI][NOD+NC0[J]] = -VAL[J];
}
if (same_string(TYP[J], "CC ") || same_string(TYP[J], "CV ") ) {
if(same_string(TYP[J], "CC ") ) T[NI][NJ]=1;
if(same_string(TYP[J], "CV ") ) T[NI][NJ]=1;
T[NI][NOD+NBR+NC0[J]] = -VAL[J];
}
}
}
/*
(3)调用GAUSS求解方程
*/
void gauss()
{
double T1, C;
for(int K = 1; K <=M; ++K) {
int I, L;
/*
这里原FORTRAN 漏掉了一些语句: I, L 没有赋值
*/
if(fabs(T[I][K])> fabs(T[L][K])) L=I;
if(fabs(T[L][K]) <1.0e-30) {
printf( "电路无唯一解 ");
exit(1);
}
if(L!=K) {
for(int J = K; J <=N; ++J) {
T1 = T[K][J];
T[K][J] = T[L][J];
T[L][J] = T1;
}
}
for (int I = K; I <=M;++I) {
C = T[I][K];
if(C!=0.0) {
for(int J = K; J <=N; ++J) {
T[I][J] = T[I][J]/C;
if(I> K) T[I][J] = T[I][J]-T[K][J];
}
}
}
}
for(int I = M-1; I> =1; --I)
for (int J = M; J> =I+1; --J)
T[I][N] = T[I][N]-T[I][J]*T[J][N];
}
/*
输出子程序
*/
void output()
{
double P;
printf( " 节点电压 \n ");
for (int J=1; J <=NOD; ++J) {
printf( "\t %d:= %gV\n ", J, T[J][N]);
}
printf( "\t 支路电压 \t 支路电流 t\支路功率\n ");
for (int J=1; J <=NBR; ++J) {
P = T[NOD+J][N]*T[NOD+NBR+J][N];
printf( "\t%d: U=%g\t I=%g\t P=%g\n ", J, T[NOD+J][N], T[NOD+NBR+J][N], P);
}
}