Quiet@Remove["Global`*"];
curl[vect_]:=Curl[vect, {r,θ,z}, "Cylindrical"];
Quiet@<<PhysicalConstants`

ϵ0=(N@VacuumPermittivity)[[1]];
μ0=(N@VacuumPermeability)[[1]];

Jvect={0,0,Jz[r]};
Hvect={0,Hfi[r],0};
rce1=Thread[curl[Jvect]==-I*ω*γ*μ*Hvect][[2]];
rce2=Thread[curl[Hvect]==Jvect][[3]];
sol=DSolve[{rce1,rce2},{Jz[r],Hfi[r]},r][[1]];
sol2=Simplify[sol,{γ>0,μ>0,ω>0}]/.C[2]->0;
jecko=Jz[r]/.sol2;
dosC1=proud==Integrate[jecko*2*Pi*r,{r,0,R}]//Solve[#,C[1]]&;
Jfinal=jecko/.dosC1[[1]];

ClearAll[Qv,Pz];
dos={γ->57.*(10^6),μ->μ0,ω->2.*Pi*f,R->0.01,delka->1,proud->1};
Qv[f_]=(1/γ)Abs[(Jfinal^2)]/.dos;
Pz[f_]:=NIntegrate[Evaluate[2*Pi*r*delka*Qv[f]/.dos],Evaluate[{r,0,R/.dos}]];
R1[f_]:=((Pz[f])/(proud^2))/.dos;
dataR={#,R1[#]}&/@Range[50,500,10];
ListPlot[dataR];

ClearAll[a];
a[f_]=(Sqrt[(2/(ω*γ*μ))])/.dos;
dat1=dataR/.{f_,R1_}:>{f,(R1/(delka/(γ*Pi*(R^2))))/.dos};
dat2=dataR/.{f_,R1_}:>{f,(R1/(delka/(γ*Pi*((R^2)-(((R-a[f]))^2)))))/.dos};
ListPlot[{dat1,dat2}]
ListPlot[dat2,PlotRange->All]

a[400]

Plot[Abs[Jfinal/.dos/.f->500],{r,0,0.01},AxesOrigin->{0,0}]

dos={γ->57.*(10^6),μ->μ0,ω->2.*Pi*f,f->50,U1->230};
a=(Sqrt[(2/(ω*γ*μ))])//.dos;
cm=0.01;
D2=20cm;
mezera=2cm;
D1=D2+2*mezera;
d1=D1+a;
d2=D2-a;
el=1;
R1=1.2*Pi*((d1*(N1^2))/(el*a*γ))/.dos;
R2=Pi*(d2/(el*a*γ))/.dos;
S1=(Pi/4)*(d1^2);
S2=(Pi/4)*(d2^2);
nagaoka=0.95;
L1=μ0*(N1^2)*(S1/el)*nagaoka;
L2=μ0*(1^2)*(S2/el)*nagaoka;
nagaokaM=0.95;
M=μ0*1*N1*(S2/el)*nagaokaM;

rcePrimar=U1==R1*i1+I*ω*L1*i1+I*ω*M*i2//.dos;
rceSekundar=0==R2*i2+I*ω*L2*i2+I*ω*M*i1//.dos;
{rcePrimar,rceSekundar};

cCu=385;
ro=8960;
V=(Pi/4)*(d1^2)*el;
T1=20;
T2=1000;
casOhrevu=15*60;
P=V*ro*cCu*(((T2-T1))/casOhrevu);
rceVykon=P==R2*Abs[(i2^2)]

resi1=Solve[rcePrimar,i1][[1]];rrr=rceVykon/.Solve[rceSekundar/.resi1,i2][[1]]
FindRoot[rrr,{N1,50}]