%Geometry.m
%%%%%%%%%%%
%Defines the geometry and the mesh of the FEMM simulation 

%Create a new magnetic problem
create(0); 
mi_probdef(SlipFreq,'meters','planar',1e-8,ActiveL,30);

%Add the three phase circuits
mi_addcircprop('A',ia,1);
mi_addcircprop('B',ib,1);
mi_addcircprop('C',ic,1);

%Drawing of the stator%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%Create origo and exterior stator nodes
mi_addnode(0,0);    %Origo
PxOut=StatorRad*sin(SlotPitch/2);
PyOut=StatorRad*cos(SlotPitch/2);
mi_addnode(PxOut,PyOut);
mi_addnode(-PxOut,PyOut);
PxIn=(RotorRad+AirGapHeight)*sin(SlotPitch/2);
PyIn=(RotorRad+AirGapHeight)*cos(SlotPitch/2);
mi_addnode(PxIn,PyIn);
mi_addnode(-PxIn,PyIn);
 
%Create the stator slot nodes
mi_clearselected();
ThetaSlotTip=(ToothOpeningWidth/2)/(RotorRad+AirGapHeight); %Angle from mid of slot to a tooth tip [rad]
PxInnerToothTip1=(RotorRad+AirGapHeight)*sin(ThetaSlotTip);
PyInnerToothTip1=(RotorRad+AirGapHeight)*cos(ThetaSlotTip);
mi_addnode(PxInnerToothTip1,PyInnerToothTip1);
mi_addnode(-PxInnerToothTip1,PyInnerToothTip1);
PxInnerToothTip2=(RotorRad+AirGapHeight+InnerToothTipHeight)*sin(ThetaSlotTip);
PyInnerToothTip2=(RotorRad+AirGapHeight+InnerToothTipHeight)*cos(ThetaSlotTip);
mi_addnode(PxInnerToothTip2,PyInnerToothTip2);
mi_addnode(-PxInnerToothTip2,PyInnerToothTip2);

ThetaSlotIn=SlotPitch/2-(ToothWidth/2)*(1/(RotorRad+AirGapHeight+OuterToothTipHeight)); %Angle from mid of slot to a tooth [rad]
PxInnerSlot=(RotorRad+AirGapHeight+OuterToothTipHeight)*sin(ThetaSlotIn);
PyInnerSlot=(RotorRad+AirGapHeight+OuterToothTipHeight)*cos(ThetaSlotIn);
mi_addnode(PxInnerSlot,PyInnerSlot);
mi_addnode(-PxInnerSlot,PyInnerSlot);
ThetaSlotOut=SlotPitch/2-(ToothWidth/2)*(1/(RotorRad+AirGapHeight+OuterToothTipHeight+ToothHeight)); %Angle from mid of slot to a tooth [rad]
PxOuterSlot=(RotorRad+AirGapHeight+OuterToothTipHeight+ToothHeight)*sin(ThetaSlotOut);
PyOuterSlot=(RotorRad+AirGapHeight+OuterToothTipHeight+ToothHeight)*cos(ThetaSlotOut);
mi_addnode(PxOuterSlot,PyOuterSlot);
mi_addnode(-PxOuterSlot,PyOuterSlot);

%Draw outer periphery arc
mi_addarc(StatorRad*sin(SlotPitch/2),StatorRad*cos(SlotPitch/2),-StatorRad*sin(SlotPitch/2),StatorRad*cos(SlotPitch/2),SlotPitch*180/pi,NoOfArcSegments);

%Draw arc corresponding to inner stator back
mi_addarc(PxOuterSlot,PyOuterSlot,-PxOuterSlot,PyOuterSlot,2*ThetaSlotIn*180/pi,NoOfArcSegments);

%Draw arc corresponding to inner tooth tip
mi_addarc(PxIn,PyIn,PxInnerToothTip1,PyInnerToothTip1,(SlotPitch/2-ThetaSlotTip)*180/pi,NoOfArcSegments);
mi_addarc(-PxInnerToothTip1,PyInnerToothTip1,-PxIn,PyIn,(SlotPitch/2-ThetaSlotTip)*180/pi,NoOfArcSegments);

%Create the stator slot lines
mi_addsegment(PxInnerToothTip1,PyInnerToothTip1,PxInnerToothTip2,PyInnerToothTip2);
mi_addsegment(-PxInnerToothTip1,PyInnerToothTip1,-PxInnerToothTip2,PyInnerToothTip2);
mi_addsegment(PxInnerToothTip2,PyInnerToothTip2,PxInnerSlot,PyInnerSlot);
mi_addsegment(-PxInnerToothTip2,PyInnerToothTip2,-PxInnerSlot,PyInnerSlot);
mi_addsegment(PxInnerSlot,PyInnerSlot,PxOuterSlot,PyOuterSlot);
mi_addsegment(-PxInnerSlot,PyInnerSlot,-PxOuterSlot,PyOuterSlot);

%Draw lines corresponding to the stator slot pitch
mi_addsegment(PxIn,PyIn,PxOut,PyOut);
mi_addsegment(-PxIn,PyIn,-PxOut,PyOut);

%Draw arc corresponding to the stator slot wedge
mi_addarc(PxInnerToothTip1,PyInnerToothTip1,-PxInnerToothTip1,PyInnerToothTip1,ThetaSlotTip*180/pi,NoOfArcSegments); %Inner part of slot wedge
mi_addarc(PxInnerToothTip2,PyInnerToothTip2,-PxInnerToothTip2,PyInnerToothTip2,2*atan(PxInnerToothTip2/PyInnerToothTip2)*180/pi,NoOfArcSegments); %Outer part of slot wedge

%Drawing of the rotor%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%Create the shaft nodes and corresponding arc
PxShaft=ShaftRad*sin(PolePitch/2);
PyShaft=ShaftRad*cos(PolePitch/2);
mi_addnode(PxShaft,PyShaft); 
mi_addnode(-PxShaft,PyShaft); 
mi_addarc(PxShaft,PyShaft,-PxShaft,PyShaft,PolePitch*180/pi,NoOfArcSegments); 

%Create the rotor slot nodes and inner rotor slot arc
PxOuterRotorSlot=OuterRotorSlotWidth/2;
PyOuterRotorSlot=RotorRad-RotorSlotDepth;
mi_addnode(PxOuterRotorSlot,PyOuterRotorSlot); 
mi_addnode(-PxOuterRotorSlot,PyOuterRotorSlot); 
PxOuterRotorSlotOpening=RotorRad*sin(RotorSlotOpeningWidth/(2*RotorRad));
PyOuterRotorSlotOpening=RotorRad*cos(RotorSlotOpeningWidth/(2*RotorRad));
mi_addnode(PxOuterRotorSlotOpening,PyOuterRotorSlotOpening); 
mi_addnode(-PxOuterRotorSlotOpening,PyOuterRotorSlotOpening); 
PxInnerRotorSlot=InnerRotorSlotWidth/2;
PyInnerRotorSlot=RotorRad-RotorSlotDepth-RotorSlotHeight;
mi_addnode(PxInnerRotorSlot,PyInnerRotorSlot); 
mi_addnode(-PxInnerRotorSlot,PyInnerRotorSlot); 
ThetaInnerRotorSlot=2*atan(PxInnerRotorSlot/PyInnerRotorSlot);
mi_addarc(PxInnerRotorSlot,PyInnerRotorSlot,-PxInnerRotorSlot,PyInnerRotorSlot,ThetaInnerRotorSlot*180/pi,NoOfArcSegments); 

%Create the rotor slot lines
mi_addsegment(PxInnerRotorSlot,PyInnerRotorSlot,PxOuterRotorSlot,PyOuterRotorSlot);
mi_addsegment(-PxInnerRotorSlot,PyInnerRotorSlot,-PxOuterRotorSlot,PyOuterRotorSlot);
mi_addsegment(PxOuterRotorSlot,PyOuterRotorSlot,PxOuterRotorSlotOpening,PyOuterRotorSlotOpening);
mi_addsegment(-PxOuterRotorSlot,PyOuterRotorSlot,-PxOuterRotorSlotOpening,PyOuterRotorSlotOpening);

%Create the outer rotor periphery nodes and corresponding arcs
PxOuterRotorMidTooth=RotorRad*sin(RotorSlotPitch/2);
PyOuterRotorMidTooth=RotorRad*cos(RotorSlotPitch/2);
mi_addnode(PxOuterRotorMidTooth,PyOuterRotorMidTooth); 
mi_addnode(-PxOuterRotorMidTooth,PyOuterRotorMidTooth); 
mi_addarc(PxOuterRotorMidTooth,PyOuterRotorMidTooth,PxOuterRotorSlotOpening,PyOuterRotorSlotOpening,0.5*(RotorSlotPitch-RotorSlotOpeningWidth/RotorRad)*180/pi,NoOfArcSegments); 
mi_addarc(-PxOuterRotorSlotOpening,PyOuterRotorSlotOpening,-PxOuterRotorMidTooth,PyOuterRotorMidTooth,0.5*(RotorSlotPitch-RotorSlotOpeningWidth/RotorRad)*180/pi,NoOfArcSegments); 
mi_addarc(-PxOuterRotorSlotOpening,PyOuterRotorSlotOpening,PxOuterRotorSlotOpening,PyOuterRotorSlotOpening,(RotorSlotOpeningWidth/RotorRad)*180/pi,NoOfArcSegments); 

%Grouping of nodes, arcs and lines and define mesh densities%%%%%%%%%%%%%%%
%Group 1: Stator
%Group 2: Shaft
%Group 3: Rotor slots

%Group stator exterior nodes
mi_clearselected();
mi_selectnode(PxOut,PyOut);
mi_selectnode(-PxOut,PyOut);
mi_selectnode(PxIn,PyIn);
mi_selectnode(-PxIn,PyIn);
mi_selectnode(PxInnerToothTip1,PyInnerToothTip1);
mi_selectnode(-PxInnerToothTip1,PyInnerToothTip1);
mi_setnodeprop('StatorExt',1); %Name the group 'StatorExt' and let it belong to group number 1

%Group the stator slot nodes
mi_clearselected();
mi_selectnode(PxOuterSlot,PyOuterSlot);
mi_selectnode(-PxOuterSlot,PyOuterSlot);
mi_selectnode(PxInnerToothTip2,PyInnerToothTip2);
mi_selectnode(-PxInnerToothTip2,PyInnerToothTip2);
mi_selectnode(PxInnerSlot,PyInnerSlot);
mi_selectnode(-PxInnerSlot,PyInnerSlot);
mi_setnodeprop('Slot',1); %Name the group 'Slot' and let it belong to group number 1

%Set properties for the arcs corresponding to the stator back
mi_clearselected();
mi_addboundprop('AZero',0, 0, 0, 0, 0, 0, 0, 0, 0); %Define the boundary condition Az=0
mi_selectarcsegment(0,StatorRad);
mi_setarcsegmentprop((MeshStatorLam/StatorRad)*180/pi,'AZero',0,1);
mi_clearselected();
mi_selectarcsegment(0,StatorRad-StatorBackThickness);
mi_setarcsegmentprop((MeshStatorLam/StatorRad)*180/pi,'UnSpec',0,1);   %UnSpec is a not defined boundary

%Set properties for the lines corresponding to the stator teeth
mi_clearselected();
mi_selectsegment((StatorRad-0.01*ToothHeight)*sin(0.5*SlotPitch),(StatorRad-0.01*ToothHeight)*cos(0.5*SlotPitch));
mi_selectsegment(-(StatorRad-0.01*ToothHeight)*sin(0.5*SlotPitch),(StatorRad-0.01*ToothHeight)*cos(0.5*SlotPitch));
mi_selectsegment((StatorRad-StatorBackThickness-0.01*ToothHeight)*sin(1*ThetaSlotIn),(StatorRad-StatorBackThickness-0.01*ToothHeight)*cos(1*ThetaSlotIn));
mi_selectsegment(-(StatorRad-StatorBackThickness-0.01*ToothHeight)*sin(1*ThetaSlotIn),(StatorRad-StatorBackThickness-0.01*ToothHeight)*cos(1*ThetaSlotIn));
mi_setsegmentprop('UnSpec',MeshStatorLam,0,0,1);    %UnSpec is a not defined boundary

%Set properties for the lines corresponding to the inner parts of the stator tooth tips and the outer part of the stator slot wedge
mi_clearselected();
PxToothTipLine=0.5*(PxInnerToothTip2+PxInnerSlot);
PyToothTipLine=PyInnerToothTip2+0.5*(PxInnerSlot-PxInnerToothTip2)*(PyInnerSlot-PyInnerToothTip2)/((PxInnerSlot-PxInnerToothTip2));
mi_selectsegment(PxToothTipLine,PyToothTipLine);
mi_selectsegment(-PxToothTipLine,PyToothTipLine);
mi_setsegmentprop('UnSpec',MeshStatorLam,0,0,1);    %UnSpec is a not defined boundary
mi_clearselected();
mi_selectarcsegment(0,RotorRad+AirGapHeight+InnerToothTipHeight);
mi_setarcsegmentprop((MeshStatorLam/(RotorRad+AirGapHeight+InnerToothTipHeight))*180/pi,'UnSpec',0,1);   %UnSpec is a not defined boundary
mi_clearselected();
mi_selectsegment((RotorRad+AirGapHeight+0.5*InnerToothTipHeight)*sin(ThetaSlotTip),(RotorRad+AirGapHeight+0.5*InnerToothTipHeight)*cos(ThetaSlotTip));
mi_selectsegment(-(RotorRad+AirGapHeight+0.5*InnerToothTipHeight)*sin(ThetaSlotTip),(RotorRad+AirGapHeight+0.5*InnerToothTipHeight)*cos(ThetaSlotTip));
mi_setsegmentprop('UnSpec',MeshStatorLam,0,0,1);    %UnSpec is a not defined boundary

%Set properties for the arcs corresponding to the outer parts of the air gap
mi_clearselected();
mi_selectarcsegment(0,RotorRad+AirGapHeight);
mi_selectarcsegment((RotorRad+AirGapHeight)*sin(SlotPitch/4+ThetaSlotTip/2),(RotorRad+AirGapHeight)*cos(SlotPitch/4+ThetaSlotTip/2));
mi_selectarcsegment(-(RotorRad+AirGapHeight)*sin(SlotPitch/4+ThetaSlotTip/2),(RotorRad+AirGapHeight)*cos(SlotPitch/4+ThetaSlotTip/2));
mi_setarcsegmentprop((MeshAirGap/(RotorRad+AirGapHeight))*180/pi,'UnSpec',0,1);   %UnSpec is a not defined boundary
 
%Group the rotor slot nodes
mi_clearselected();
mi_selectnode(PxOuterRotorSlot,PyOuterRotorSlot);
mi_selectnode(-PxOuterRotorSlot,PyOuterRotorSlot);
mi_selectnode(PxInnerRotorSlot,PyInnerRotorSlot);
mi_selectnode(-PxInnerRotorSlot,PyInnerRotorSlot);
mi_selectnode(PxOuterRotorSlotOpening,PyOuterRotorSlotOpening);
mi_selectnode(-PxOuterRotorSlotOpening,PyOuterRotorSlotOpening);
mi_setnodeprop('RotorSlot',3); %Name the group 'RotorSlot' and let it belong to group number 3
 
%Set properties for the arcs corresponding to the inner part of the air gap belonging to the rotor 
mi_clearselected();
mi_selectarcsegment(0,RotorRad);
mi_selectarcsegment(RotorRad*sin(RotorSlotPitch/4+RotorSlotOpeningWidth/(4*RotorRad)),RotorRad*cos(RotorSlotPitch/4+RotorSlotOpeningWidth/(4*RotorRad))); 
mi_selectarcsegment(-RotorRad*sin(RotorSlotPitch/4+RotorSlotOpeningWidth/(4*RotorRad)),RotorRad*cos(RotorSlotPitch/4+RotorSlotOpeningWidth/(4*RotorRad))); 
mi_setarcsegmentprop((MeshAirGap/RotorRad)*180/pi,'UnSpec',0,3);   %UnSpec is a not defined boundary

%Set properties for the inner arc of the rotor slot
mi_clearselected();
mi_selectarcsegment(0,RotorRad-RotorSlotDepth-RotorSlotHeight);
mi_setarcsegmentprop((MeshRotor/(RotorRad-RotorSlotDepth-RotorSlotHeight))*180/pi,'UnSpec',0,3);   %UnSpec is a not defined boundary

%Set properties for the rotor slot lines
mi_clearselected();
mi_selectsegment((PxOuterRotorSlot+PxInnerRotorSlot)/2,(PyOuterRotorSlot+PyInnerRotorSlot)/2);
mi_selectsegment(-(PxOuterRotorSlot+PxInnerRotorSlot)/2,(PyOuterRotorSlot+PyInnerRotorSlot)/2);
mi_selectsegment((PxOuterRotorSlot+PxOuterRotorSlotOpening)/2,(PyOuterRotorSlot+PyOuterRotorSlotOpening)/2);
mi_selectsegment(-(PxOuterRotorSlot+PxOuterRotorSlotOpening)/2,(PyOuterRotorSlot+PyOuterRotorSlotOpening)/2);
mi_setsegmentprop('UnSpec',MeshRotor,0,0,3);    %UnSpec is a not defined boundary

%Set properties for the arc corresponding to the shaft
mi_clearselected();
mi_selectarcsegment(0,ShaftRad);
mi_setarcsegmentprop((MeshShaft/ShaftRad)*180/pi,'UnSpec',0,2);   %UnSpec is a not defined boundary

%Rotational copy of stator and rotor%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
mi_clearselected();
mi_seteditmode('group');
mi_selectgroup(1);  %Select the stator
mi_copyrotate(0,0,SlotPitch*180/pi,Qs-1);

%Rotational copy of the shaft
mi_clearselected();
mi_seteditmode('group');
mi_selectgroup(2);  %Select the rotor
mi_copyrotate(0,0,PolePitch*180/pi,p-1);

%Rotational copy of the rotor (except the shaft)
mi_clearselected();
mi_seteditmode('group');
mi_selectgroup(3);  %Select the rotor
mi_copyrotate(0,0,RotorSlotPitch*180/pi,Qr-1);