function critical_current_density_nb3sn(c_pt,T,B,St,Tc0,Bc) result(critical_current_density) bind(C)
type(c_ptr), value :: c_pt
real(dp), value :: T,B,St,Tc0,Bc
real(dp) :: critical_current_density
real(dp) :: ee,s,tt,bb,h,fp,Blim,Tlim,Blow,Tlow
type(nb3sn_cfg_t), pointer :: nb3sn_ptr
call c_f_pointer(c_pt,nb3sn_ptr)
! L. Bottura, B. Bordini, Jc(B,T,e) Parameterization for the ITER
! Nb3Sn Production, IEEE Trans. Appl. Sup., 19(2), 1477-1480, 2009
Blow=1.0e-3_dp
Tlow=0.0_dp
Blim=max(abs(B),Blow)
Tlim=max(T,Tlow)
ee = St-nb3sn_ptr%emax
s = sNb3Sn(ee,nb3sn_ptr%Ca1,nb3sn_ptr%Ca2,nb3sn_ptr%e0a)
tt = Tlim/Tc0
if(tt>1.0_dp) then
critical_current_density = 0.0_dp
return
endif
bb = Blim/Bc
if(bb>1.0_dp) then
critical_current_density = 0.0_dp
return
endif
h = hNb3Sn(tt,nb3sn_ptr%nu)
fp = fpNb3Sn(bb,nb3sn_ptr%p,nb3sn_ptr%q)
critical_current_density = nb3sn_ptr%C0/Blim * s * h * fp
end function critical_current_density_nb3sn