/out,m_ind,out
/com,m_ind.mac
!  The rotor is rotated in a CCW direction and the inductance maxtrix is calculated
!  for each position of the rotor.
!
!    arg1 = starting mechanical angle (defaults to zero)
1
!    arg2 = ending mechanical angle (defaults to the edge of the model)
1
!    arg3 = current to be applied to the coils
1
!    arg4 = equivalent current for the magnet
!
!    arg5 = total number of solutions. The total number
!           of steps also includes the starting and the ending angles. 
!           (defaults to 1, only one solution, at the starting angle).
!           number of intermediate solution must be less than 254
!
!    arg6 = name of the plot file (Must be enclosed in single quotes)
!           The extention of the file  .f33  The name defaults to
!           "mname".f33
!           This is also the name of the ASCII file containing the
!           inductances for each angle, "arg6".lma or "mname".lma
!
!    arg7 = 0,1 generate BSUM plot for each operating solution which
!           is the first step for the inductance calcalation
!         = 2 generate  B vector plot for each solution
!         = -1 No plots are generated
!           The plots are generated according to the settings
!           in the user supplied macro  'winset.mac'  This macro
!           is generated by obtaining the view of the core as
!           as needed by the graphics commands and then using
!           /gsave,winset,mac to generate the macro winset.mac         
!
!    arg9 = 1, user is prompted for the first five arguements
!
!  output:
!     array of inductances: M_L_MAT(i)
!         Ith      Inductance (H)
!          1           Laa
!          2           Lbb
!          3           Lcc
!          4           Lff
!          5           Lab
!          6           Lac
!          7           Laf
!          8           Lbc
!          9           Lbf
!         10           Lcf
!         11           Angle of rotor
!
!    Format of the ASCII file containing the Inductance array
!
!    This file is in three parts. All numbers are in F10.x format
!
!   Part I consists of a header record of one line
!        Column         Definition
!          1          Number of solutions
!          2          Current for the magnet
!          3          Current for the coils
!       
!   Part II consists of the diagonal terms, Laa,Lbb,Lcc  (mH)
!        Column     Definition
!          1       Rotor angle
!          2           Laa
!          3           Lbb
!          4           Lcc
!
!   Part III consists of the off-diagonal terms (mH)
!          1       Rotor angle
!          2           Lff
!          3           Lab
!          4           Lac
!          5           Laf
!          6           Lbc
!          7           Lbf
!          8           Lcf
!     
!    macros called:
!     MVROT3D   moves the rotor and connects the rotor to the stator
!               by constraint equations (for 3d)
!     MVROTOR   moves the rotor and connects the rotor to the stator
!               by constraint equations (for 2d)
!     LD_COIL   applies the current to the coils
!     MACHMATX  Compute the inductance matrix

!
!
!   Assumptions:
!    The winding was constructed using the b_wndsc.mac
!    Phases are equally spaced in terms of electrical angle
!    (three phases are 120 degrees apart)
!
!
!   parameter file of all data:   m_ind.par
!
bmx_plt=2.6   !  maximum BSUM value for BSUM plot
!
!
*if,arg9,eq,1,then
 *ask,_arg_1,Starting Mechanical Angle,0
 *ask,_arg_2,Ending Mechanical Angle,0
 *ask,_arg_3,Coil Current (A),10
 *ask,_arg_4,Magnet Current,10
 *ask,_arg_5,Number of solutions,10
 arg1=_arg_1
 arg2=_arg_2
 arg3=_arg_3
 arg4=_arg_4
 arg5=_arg_5
 arg9=0
*endif

/nerr,0,1e5
esel,none
cmsel,,%s_iron_c%
*get,_ecnt,elem,,count
_err=0
*if,_ecnt,le,0,then
 cmsel,,s_iron
 *get,_ecnt,elem,,count
 *if,_ecnt,le,0,then
  /nerr
  /out
  *msg,error
  The element component for the iron is not defined as&
  S_IRON or the parameter S_IRON_C has not been defined&
  No action
  /out,m_ind,out,,append
  _err=1
 *else
   s_iron_c='s_iron'
 *endif
*endif

/nerr,0,1e5
*if,_err,eq,1,:end

dsys

*get,mcpu1,active,,time,cpu
*if,arg1,le,0,then
   arg1=0
*endif
_stang=arg1

_endang=arg2
*if,arg2,le,0,then
  esel,all
  esel,u,ename,,36
  nsle
  csys,1
  *get,_endang,node,,mxloc,y
*endif

_arg6=arg6
*get,_targ6,parm,_arg6,type
*if,_targ6,ne,3,then
  _arg6=mname
*endif

*if,arg3,le,0,then
   arg3=1
*endif
_pkamp=arg3

/com, apply the force boundary conditions-name of iron
/com, is obtained from b_wndsc.mac
_nomsg=1 !  turn off the gui summary from macros        
_mg1=1   !  turn off the gui summary from macros        


/psf,mxwf,1
nsel,all
esel,all
/pnum,node
nplo
*if,arg5,le,0,then
  _numcp=1
*else
  _numcp=arg5
*endif
*if,_numcp,gt,256,then
  _numcp=256
*endif

cmsel,,%s_iron_c%
nsle
nsel,r,ext
csys,1
*get,_mxang,node,,mxloc,y
*get,_mnang,node,,mnloc,y

_fulmod=0
*if,_mxang-_mnang,gt,180.1,then
  _fulmod=1
*endif

*get,_dimn,active,,solu,dimn        ! 1=axisym, 2=planar, 3=3d 
*if,arg7,ne,-1,then
 *if,_dimn,eq,3,then
  winset
  _plt_cur=0
 *else
  _plt_cur=1   !  supimpose the currents
  mvrotor,arg2
  /out,m_ind,out,,append
  /sho,xxx
  /auto
  /plopt,info,on
  alls
  eplo
  /sho,term
  eplo
  /user
  /plopt,wins,off
 *endif
 /sho,%_arg6%,f33
*endif

cmsel,,%s_iron_c%
nsel,,ext
*get,mxn_s,node,,num,max

!   electrical angle for each solution
*set,_rphst
*dim,_rphst,,_numcp
*if,_numcp,gt,1,then
 _anginc=(_endang-_stang)/(_numcp-1)
*else
 _anginc=0
*endif
_c_ang=_stang-_anginc
*set,_timcyc
*dim,_timcyc,,_numcp

/com, compute the periodic factor for the axial
esel,,ename,,36
nsle
*get,_zcmx,node,,mxloc,z
*get,_zcmn,node,,mnloc,z
cmsel,,%s_iron_c%
nsle
*get,_zsmx,node,,mxloc,z
*get,_zsmn,node,,mnloc,z

*get,_dimn,active,,solu,dimn        ! 1=axisym, 2=planar, 3=3d 
*if,abs(_zsmx-_zsmn),lt,.001,then
 z_factor=1
 *if,_dimn,ne,3,then
  z_factor=stkthk*ggeom
 *endif
*else
 z_factor=nint((_zcmx-_zcmn)/(_zsmx-_zsmn))
 *if,z_factor,gt,2,then
   z_factor=2
 *endif
*endif
*set,m_l_ind
*dim,m_l_ind,,_numcp,11

*do,_ix1,1,_numcp
 /gopr
 *get,_t1cyc1,active,,time,cpu
 *msg,info,_ix1
 Cycle in _ix1:  %i
 _c_ang=_c_ang+_anginc

 _rphst(_ix1)=_c_ang

 *if,_fulmod,eq,1,then
  _anginp=_anginc
 *else
  _anginp=_c_ang
 *endif
 _nomsg=1
 machmatx,_anginp,0,arg3,arg4
 /out,m_ind,out,,append
 *msg,info,_anginp,arg3,arg4
 _anginp,arg3,arg4:  %g %g %g
 m_l_ind(_ix1,11)=_c_ang
 *do,_ix2,1,10
  m_l_ind(_ix1,_ix2)=l_matx(_ix2)
 *enddo

 *get,_t2cyc2,active,,time,cpu
 _timcyc(_ix1)=_t2cyc2-_t1cyc1
 parsav,all,m_ind,par
*enddo

*get,mcpu2,active,,time,cpu
mdcpu=mcpu2-mcpu1

*get,_dimn,active,,solu,dimn        ! 1=axisym, 2=planar, 3=3d 
*if,_dimn,eq,3,then
 ld_coil,_c_eang,0   !  reset the coil real constant sets
*endif

/out,m_ind,out,,append
/nerr,0,1e5
/nopr
/out,m_ind,sum
*msg,info
_________SUMMARY OF SOLUTIONS____________
*msg,info,mname
Model ID:_________________________ %c
*msg,info,_arg6
Plot File name:___________________ %c .f33
*msg,info,_arg6
Table of Inductances File name:___ %c .lma (ASCII file)
*msg,info,arg1
Starting mechanical angle(D):_____ %g
*msg,info,arg2
Ending mechanical angle(D):_______ %g
*msg,info,arg3
Peak Current applied (A):_________ %g
*msg,info,arg4
Magnet Current(A):________________ %g
*msg,info,n_wnd_f
Winding File:_____________________ %c
*msg,info,arg5
Total number of solutions:________ %g
*msg,info,s_iron_c
Name of stator component:_________ %c
*if,_dimn,ne,3,then
 *msg,info,stkthk
 Stack thickness:__________________ %g  (2D model)
*else
 *msg,info,z_factor
 Axial periodic factor:____________ %g (3D model)
*endif
*msg,info,a_factor
 Angular periodicity factor:_______ %g

*vwrit,mdcpu
(' Total CPU time (sec):_____________ ',f6.0)
*msg,info


*msg,info
 ______Mechanical_______Diagonal (Lii)_______
*msg,info
 _______Angle (D)_________ (mH)______________
*msg,info
 ______________Laa_____Lbb_____Lcc_____Lff___
!23456789012345678901234567890!23456789012345678901234567890!2345678901234567890
!        1!        1!        1!        1!        1!        1!        1!        1
!           xxxxxx.yyzzzzzz.yyzzzzzz.yy

*vwrite,m_l_ind(1,11),m_l_ind(1,1),m_l_ind(1,2),m_l_ind(1,3),m_l_ind(1,4)
(f12.2,4f8.2)
*msg,info

*msg,info

!23456789012345678901234567890!23456789012345678901234567890!2345678901234567890
!        1!        1!        1!        1!        1!        1!        1!        1
!           xxxxxx.yyzzzzzz.yyzzzzzz.yy
*msg,info
 ______Mechanical______________OFF_Diagonal (Lij)__________________
*msg,info
 _______Angle (D)___________________ (mH)__________________________
*msg,info
 ______________Lab_____Lac_____Laf_____Lbc_____Lbf_____Lcf_________
*abbr,_vwr1,*vwrite,m_l_ind(1,11),m_l_ind(1,5),m_l_ind(1,6)
*vlen,18
_vwr1,m_l_ind(1,7),m_l_ind(1,8),m_l_ind(1,9),m_l_ind(1,10)
(f12.2,7f8.2)

*msg,info
___________________________________________________________________
*msg,info

*msg,info
Inductance is for a 360 degree model and&
for the total thickness of the stack
*msg,info

*msg,info
The tabular results in the ASCII file are in the&
same order as listed here

/out,m_ind,out,,append
/gopr
*set,_aval
*dim,_aval,,1,3
_aval(1,1)=arg5
_aval(1,2)=arg3
_aval(1,3)=arg4
m_l_ind(1,11)=m_l_ind(1,11)+1e-5
/nopr
/com, write out the ASCII file to *.lma
/out,%_arg6%,lma
*vwrit,_aval(1,1),_aval(1,2),_aval(1,3)
(f10.0,2f10.5)
*vwrite,m_l_ind(1,11),m_l_ind(1,1),m_l_ind(1,2),m_l_ind(1,3)
(f10.3,3f10.4)
*abbr,_vwr1,*vwrit,m_l_ind(1,11),m_l_ind(1,4),m_l_ind(1,5)
_vwr1,m_l_ind(1,6),m_l_ind(1,7),m_l_ind(1,8),m_l_ind(1,9),m_l_ind(1,10)
(f10.3,7f10.4)

/out,machmatx,out,,append

*uili,m_ind,sum
*abbr,_vwr1
:end

/sho,xxx
/sho,term
/out

/gopr
_mg1=0
_nomsg=0