>>740
#include <stdio.h>
#include <math.h>
double f(double B, double k2, double k1, double T){
double b;
b=tan(sqrt(k1*k1-B*B)*T)-sqrt((B*B-k2*k2)/(k1*k1-B*B));
return b;}
int main(void){
int p;
double k1,k2,pai=3.141592,r0=0.03;
double b,T,a,S,R,d,ep=1.56,Bet;//R(Xn-1) S(Xn) Bet(Xn+1) ep(ε)
int N=500,i,j;
for(j=1;j<N+1;j++){
T=0.01*sin(pai*j*0.01)+0.011;
k2=2*pai/r0;
k1=ep*k2;
for(b=k2;b<k1;b++){
a=f(b,k2,k1,T);
if(a<0){R=b; break;}
else { S=b; }}
for(p=0;p<10;p++){
d=f(S,k2,k1,T)-f(R,k2,k1,T);
if(d==0){break;}
Bet=(R*f(S,k2,k1,T)-S*f(R,k2,k1,T))/(f(S,k2,k1,T)-f(R,k2,k1,T));
R=S,S=Bet;}
// 〜中略}
return 0;}

bとかサブルーチンに関わる部分はこのような感じです
よろしくお願いします。