In this research, the electrothermal characteristics of anode have been analytically investigated by considering two subzones including the quasi-neutral and non-neutral plasma flow near the anode surface. A correlation between current density and anode fall voltage has been derived by using the continuity and momentum equations for electron in the non-neutral plasma subzone. Also, two relations have been deduced from the magnetic field induction equation and Ohm’s law to determine the current densities and voltage drops in the quasi-neutral plasma subzone. The suggested solution algorithm has been based on the satisfaction of current continuity condition to be converged by a trial and error method. The Princeton benchmark thruster has been investigated for the discharge current of 8 kA with mass flow rates of 4 and 24 gr/s to evaluate the developed algorithm. The obtained analytical results show that under operating conditions with mass flow rates of 4 and 24 g/s, the heat flux values are respectively in the range of 480 to1350 and 400 to1490 W/cm2, the current densities vary from 24 to 90 and from 33 to 140 A/cm2, and also the anode falls change from 5.8 to 11 and 4.9 to 6 V. The analytical algorithm has predicted the heat fluxes, current densities, and anode falls along the anode face in good agreement with measurements. It has been seen that the heat flux and current density around the anode mid-lip, where the anode fall has its minimum value, have been maximized. By using the second thermodynamic law, it has been shown that the probability density function of effective velocity has a Gaussian distribution. Therefore, under this circumstance, entropy has been maximized, and consequently the heat flux and current density profiles have been taken the summit form.