📄 femmviewdoc.cpp
字号:
v=v+i+1;
i=k;
}
k=strlen(v);
if(k>0) for(i=k-1;i>=0;i--){
if(v[i]=='\"'){
v[i]=0;
i=-1;
}
}
MProp.BlockName=v;
q[0]=NULL;
}
if( _strnicmp(q,"<mu_x>",6)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.mu_x);
q[0]=NULL;
}
if( _strnicmp(q,"<mu_y>",6)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.mu_y);
q[0]=NULL;
}
if( _strnicmp(q,"<H_c>",5)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.H_c);
q[0]=NULL;
}
if( _strnicmp(q,"<J_re>",6)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Jr);
q[0]=NULL;
}
if( _strnicmp(q,"<J_im>",6)==0){
v=StripKey(s);
if (Frequency!=0) sscanf(v,"%lf",&MProp.Ji);
q[0]=NULL;
}
if( _strnicmp(q,"<sigma>",7)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Cduct);
q[0]=NULL;
}
if( _strnicmp(q,"<phi_h>",7)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Theta_hn);
q[0]=NULL;
}
if( _strnicmp(q,"<phi_hx>",8)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Theta_hx);
q[0]=NULL;
}
if( _strnicmp(q,"<phi_hy>",8)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Theta_hy);
q[0]=NULL;
}
if( _strnicmp(q,"<d_lam>",7)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.Lam_d);
q[0]=NULL;
}
if( _strnicmp(q,"<LamFill>",8)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.LamFill);
q[0]=NULL;
}
if( _strnicmp(q,"<LamType>",9)==0){
v=StripKey(s);
sscanf(v,"%i",&MProp.LamType);
q[0]=NULL;
}
if( _strnicmp(q,"<NStrands>",10)==0){
v=StripKey(s);
sscanf(v,"%i",&MProp.NStrands);
q[0]=NULL;
}
if( _strnicmp(q,"<WireD>",7)==0){
v=StripKey(s);
sscanf(v,"%lf",&MProp.WireD);
q[0]=NULL;
}
if( _strnicmp(q,"<BHPoints>",10)==0){
v=StripKey(s);
sscanf(v,"%i",&MProp.BHpoints);
if (MProp.BHpoints>0)
{
MProp.Hdata=(CComplex *)calloc(MProp.BHpoints,sizeof(CComplex));
MProp.Bdata=(double *)calloc(MProp.BHpoints,sizeof(double));
for(j=0;j<MProp.BHpoints;j++){
fgets(s,1024,fp);
MProp.Hdata[j]=0;
sscanf(s,"%lf %lf",&MProp.Bdata[j],&MProp.Hdata[j].re);
}
}
q[0]=NULL;
}
if( _strnicmp(q,"<endblock>",9)==0){
if (MProp.BHpoints>0) MProp.GetSlopes(Frequency*2.*PI);
blockproplist.Add(MProp);
MProp.BHpoints=0;
MProp.Bdata=NULL;
MProp.Hdata=NULL;
MProp.slope=NULL;
q[0]=NULL;
}
// Circuit Properties
if( _strnicmp(q,"<begincircuit>",14)==0){
CProp.CircName="New Circuit";
CProp.CircType=0;
CProp.Amps=0;
q[0]=NULL;
}
if( _strnicmp(q,"<circuitname>",13)==0){
v=StripKey(s);
k=strlen(v);
for(i=0;i<k;i++)
if(v[i]=='\"'){
v=v+i+1;
i=k;
}
k=strlen(v);
if(k>0) for(i=k-1;i>=0;i--){
if(v[i]=='\"'){
v[i]=0;
i=-1;
}
}
CProp.CircName=v;
q[0]=NULL;
}
if( _strnicmp(q,"<totalamps_re>",14)==0){
double inval;
v=StripKey(s);
sscanf(v,"%lf",&inval);
CProp.Amps+=inval;
q[0]=NULL;
}
if( _strnicmp(q,"<totalamps_im>",14)==0){
double inval;
v=StripKey(s);
sscanf(v,"%lf",&inval);
if (Frequency!=0) CProp.Amps+=(I*inval);
q[0]=NULL;
}
if( _strnicmp(q,"<circuittype>",13)==0){
v=StripKey(s);
sscanf(v,"%i",&CProp.CircType);
q[0]=NULL;
}
if( _strnicmp(q,"<endcircuit>",12)==0){
circproplist.Add(CProp);
q[0]=NULL;
}
// Points list;
if(_strnicmp(q,"[numpoints]",11)==0){
v=StripKey(s);
sscanf(v,"%i",&k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%lf %lf %i\n",&node.x,&node.y,&t);
node.BoundaryMarker=t-1;
nodelist.Add(node);
}
q[0]=NULL;
}
// read in segment list
if(_strnicmp(q,"[numsegments]",13)==0){
v=StripKey(s);
sscanf(v,"%i",&k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%i %i %lf %i %i\n",&segm.n0,&segm.n1,&segm.MaxSideLength,&t,&segm.Hidden);
segm.BoundaryMarker=t-1;
linelist.Add(segm);
}
q[0]=NULL;
}
// read in arc segment list
if(_strnicmp(q,"[numarcsegments]",13)==0){
v=StripKey(s);
sscanf(v,"%i",&k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%i %i %lf %lf %i %i\n",&asegm.n0,&asegm.n1,
&asegm.ArcLength,&asegm.MaxSideLength,&t,&asegm.Hidden);
asegm.BoundaryMarker=t-1;
arclist.Add(asegm);
}
q[0]=NULL;
}
// read in list of holes;
if(_strnicmp(q,"[numholes]",13)==0){
v=StripKey(s);
sscanf(v,"%i",&k);
if(k>0)
{
blk.BlockType=-1;
blk.MaxArea=0;
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%lf %lf\n",&blk.x,&blk.y);
// blocklist.Add(blk);
// don't add holes to the list
// of block labels because it messes up the
// number of block labels.
}
}
q[0]=NULL;
}
// read in regional attributes
if(_strnicmp(q,"[numblocklabels]",13)==0){
v=StripKey(s);
sscanf(v,"%i",&k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%lf %lf %i %lf %i %lf %i %i %i\n",&blk.x,&blk.y,
&blk.BlockType,&blk.MaxArea,
&blk.InCircuit,&blk.MagDir,
&blk.InGroup,&blk.Turns,
&blk.IsExternal);
if (blk.MaxArea<0) blk.MaxArea=0;
else blk.MaxArea=PI*blk.MaxArea*blk.MaxArea/4.;
blk.BlockType-=1;
blk.InCircuit-=1;
blocklist.Add(blk);
}
q[0]=NULL;
}
if(_strnicmp(q,"[solution]",10)==0){
flag=TRUE;
q[0]=NULL;
}
}
// read in meshnodes;
fscanf(fp,"%i\n",&k);
meshnode.SetSize(k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
if (Frequency!=0)
sscanf(s,"%lf %lf %lf %lf",&mnode.x,&mnode.y,
&mnode.A.re,&mnode.A.im);
else{
sscanf(s,"%lf %lf %lf",&mnode.x,&mnode.y,&mnode.A.re);
mnode.A.im=0;
}
meshnode.SetAt(i,mnode);
}
// read in elements;
fscanf(fp,"%i\n",&k);
meshelem.SetSize(k);
for(i=0;i<k;i++)
{
fgets(s,1024,fp);
sscanf(s,"%i %i %i %i",&elm.p[0],&elm.p[1],&elm.p[2],&elm.lbl);
elm.blk=blocklist[elm.lbl].BlockType;
meshelem.SetAt(i,elm);
}
// read in circuit data;
fscanf(fp,"%i\n",&k);
for(i=0;i<k;i++){
fgets(s,1024,fp);
if (Frequency==0){
sscanf(s,"%i %lf",&j,&zr);
blocklist[i].Case=j;
if (j==0) blocklist[i].dVolts=zr;
else blocklist[i].J=zr;
}
else{
sscanf(s,"%i %lf %lf",&j,&zr,&zi);
blocklist[i].Case=j;
if (j==0) blocklist[i].dVolts=zr + I*zi;
else blocklist[i].J=zr + I*zi;
}
}
fclose(fp);
// scale depth to meters for internal computations;
if(Depth==-1) Depth=1; else Depth*=LengthConv[LengthUnits];
// Find flux density in each element;
for(i=0;i<meshelem.GetSize();i++) GetElementB(meshelem[i]);
// Find extreme values of A;
A_Low=meshnode[0].A.re; A_High=meshnode[0].A.re;
for(i=1;i<meshnode.GetSize();i++)
{
if (meshnode[i].A.re>A_High) A_High=meshnode[i].A.re;
if (meshnode[i].A.re<A_Low) A_Low =meshnode[i].A.re;
if(Frequency!=0)
{
if (meshnode[i].A.im<A_Low) A_Low =meshnode[i].A.im;
if (meshnode[i].A.im>A_High) A_High=meshnode[i].A.im;
}
}
// save default values for extremes of A
A_lb=A_Low;
A_ub=A_High;
if(Frequency!=0){ // compute frequency-dependent permeabilities for linear blocks;
CComplex deg45; deg45=1+I;
CComplex K,halflag;
double ds;
double w=2.*PI*Frequency;
for(k=0;k<blockproplist.GetSize();k++){
if (blockproplist[k].LamType==0){
blockproplist[k].mu_fdx = blockproplist[k].mu_x*
exp(-I*blockproplist[k].Theta_hx*PI/180.);
blockproplist[k].mu_fdy = blockproplist[k].mu_y*
exp(-I*blockproplist[k].Theta_hy*PI/180.);
if(blockproplist[k].Lam_d!=0){
halflag=exp(-I*blockproplist[k].Theta_hx*PI/360.);
ds=sqrt(2./(0.4*PI*w*blockproplist[k].Cduct*blockproplist[k].mu_x));
K=halflag*deg45*blockproplist[k].Lam_d*0.001/(2.*ds);
blockproplist[k].mu_fdx=(blockproplist[k].mu_fdx*tanh(K)/K)*
blockproplist[k].LamFill+(1.-blockproplist[k].LamFill);
halflag=exp(-I*blockproplist[k].Theta_hy*PI/360.);
ds=sqrt(2./(0.4*PI*w*blockproplist[k].Cduct*blockproplist[k].mu_y));
K=halflag*deg45*blockproplist[k].Lam_d*0.001/(2.*ds);
blockproplist[k].mu_fdy=(blockproplist[k].mu_fdy*tanh(K)/K)*
blockproplist[k].LamFill+(1.-blockproplist[k].LamFill);
}
}
}}
// element centroids and radii;
for(i=0;i<meshelem.GetSize();i++)
{
meshelem[i].ctr=Ctr(i);
for(j=0,meshelem[i].rsqr=0;j<3;j++)
{
b=sqr(meshnode[meshelem[i].p[j]].x-meshelem[i].ctr.re)+
sqr(meshnode[meshelem[i].p[j]].y-meshelem[i].ctr.im);
if(b>meshelem[i].rsqr) meshelem[i].rsqr=b;
}
}
// build list of elements connected to each node;
// allocate connections list;
NumList=(int *)calloc(meshnode.GetSize(),sizeof(int));
ConList=(int **)calloc(meshnode.GetSize(),sizeof(int *));
// find out number of connections to each node;
for(i=0;i<meshelem.GetSize();i++)
for(j=0;j<3;j++)
NumList[meshelem[i].p[j]]++;
// allocate space for connections lists;
for(i=0;i<meshnode.GetSize();i++)
ConList[i]=(int *)calloc(NumList[i],sizeof(int));
// build list;
for(i=0;i<meshnode.GetSize();i++) NumList[i]=0;
for(i=0;i<meshelem.GetSize();i++)
for(j=0;j<3;j++){
k=meshelem[i].p[j];
ConList[k][NumList[k]]=i;
NumList[k]++;
}
// find extreme values of J;
{
CComplex Jelm[3],Aelm[3];
GetJA(0,Jelm,Aelm);
Jr_Low=fabs(Jelm[0].re);
Jr_High=Jr_Low;
Ji_Low=fabs(Jelm[0].im);
Ji_High=Ji_Low;
J_Low=abs(Jelm[0]);
J_High=J_Low;
for(i=0;i<meshelem.GetSize();i++)
{
GetJA(i,Jelm,Aelm);
for(j=0;j<3;j++){
br=fabs(Jelm[j].re);
bi=fabs(Jelm[j].im);
b=abs(Jelm[j]);
if(b>J_High) J_High=b; if(b<J_Low) J_Low=b;
if(br>Jr_High) Jr_High=br; if(br<Jr_Low) Jr_Low=br;
if(bi>Ji_High) Ji_High=bi; if(bi<Ji_Low) Ji_Low=bi;
}
}
}
// Find extreme values of B;
Br_Low=sqrt(sqr(meshelem[0].B1.re) + sqr(meshelem[0].B2.re));
Br_High=Br_Low;
Bi_Low=sqrt(sqr(meshelem[0].B1.im)+ sqr(meshelem[0].B2.im));
Bi_High=Bi_Low;
B_Low=sqrt(Br_Low*Br_Low + Bi_Low*Bi_Low);
B_High=B_Low;
for(i=0;i<meshelem.GetSize();i++)
{
GetNodalB(meshelem[i].b1,meshelem[i].b2,meshelem[i]);
for(j=0;j<3;j++){
br=sqrt(sqr(meshelem[i].b1[j].re) +
sqr(meshelem[i].b2[j].re));
bi=sqrt(sqr(meshelem[i].b1[j].im) +
sqr(meshelem[i].b2[j].im));
b=sqrt(br*br+bi*bi);
if(b>B_High) B_High=b; if(b<B_Low) B_Low=b;
if(br>Br_High) Br_High=br; if(br<Br_Low) Br_Low=br;
if(bi>Bi_High) Bi_High=bi; if(bi<Bi_Low) Bi_Low=bi;
}
}
// Find extreme values of H;
H_High=0;
for(i=0;i<meshelem.GetSize();i++)
{
if(Frequency!=0)
{
CComplex b1,b2,h1,h2,mu1,mu2;
double h;
b1=meshelem[i].B1;
b2=meshelem[i].B2;
GetMu(b1,b2,mu1,mu2,i);
h1 = b1/(mu1*muo);
h2 = b2/(mu2*muo);
h = sqrt(Re(h1*conj(h1)+h2*conj(h2)));
if (h>H_High) H_High=h;
}
else{
double b1,b2,h1,h2,mu1,mu2;
double h;
b1=Re(meshelem[i].B1);
b2=Re(meshelem[i].B2);
GetMu(b1,b2,mu1,mu2,i);
h1 = b1/(mu1*muo);
h2 = b2/(mu2*muo);
h = sqrt(h1*h1+h2*h2);
if (h>H_High) H_High=h;
}
}
// Choose bounds based on the type of contour plot
// currently in play
POSITION pos = GetFirstViewPosition();
CFemmviewView *theView=(CFemmviewView *)GetNextView(pos);
if(Frequency==0)
{
if (theView->DensityPlot==2) theView->DensityPlot=1;
if (theView->DensityPlot>1) theView->DensityPlot=0;
}
if (theView->DensityPlot==0){ B_lb=B_Low; B_ub=B_High;}
if (theView->DensityPlot==1){ B_lb=B_Low; B_ub=B_High;}
if (theView->DensityPlot==2){ B_lb=Br_Low; B_ub=Br_High;}
if (theView->DensityPlot==3){ B_lb=Bi_Low; B_ub=Bi_High;}
if (theView->DensityPlot==4){ B_lb=J_Low; B_ub=J_High;}
if (theView->DensityPlot==5){ B_lb=Jr_Low; B_ub=Jr_High;}
if (theView->DensityPlot==6){ B_lb=Ji_Low; B_ub=Ji_High;}
// compute total resulting current for circuits with an a priori defined
// voltage gradient; Need this to display circuit results & impedance.
for(i=0;i<circproplist.GetSize();i++)
{
CComplex Jelm[3],Aelm[3];
double a;
if(circproplist[i].CircType>1)
for(j=0,circproplist[i].Amps=0.;j<meshelem.GetSize();j++)
{
if(blocklist[meshelem[j].lbl].InCircuit==i)
{
GetJA(j,Jelm,Aelm);
a=ElmArea(j)*sqr(LengthConv[LengthUnits]);
for(k=0;k<3;k++) circproplist[i].Amps+=a*Jelm[k]/3;
}
}
}
// Build adjacency information for each element.
// This info is required to find points within the
// mesh quickly, instead of by a brute-force search.
BuildElemMap(&meshnode, &meshelem, NumList, ConList, (void *)this);
// Check to see if any regions are multiply defined
// (i.e. tagged by more than one block label). If so,
// display an error message and mark the problem blocks.
for(k=0,bMultiplyDefinedLabels=FALSE;k<blocklist.GetSize();k++)
{
if((i=InTriangle(blocklist[k].x,blocklist[k].y))>=0)
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -