I would like to implement a solver (in C) for the Diophantine equation $x^2+3y^2=N$ for non-negative integers $\{x,y\}$ and positive integer $N$. I have read online that one has to prime factorize N and express the primes in $\mathbb{Z}[\sqrt{-3}]$. Unfortunately I do not have sufficent mathematical expertise to follow through the heavy proofs for the last part.
Is there an easy-to-understand method of expressing primes in $\mathbb{Z}[\sqrt{-3}]$ that can also be efficiently implemented in software?
PS: My background in mathematics consists of engineering calculus and basic number theory up to modular arithmetic.
I guess in your case it might be useful to learn about Eisenstein integers. They are given by $$ z= a+ b \omega$$ with $a$ and $b$ integers and $\omega = \exp(2\pi i/3) = (i\sqrt{3} -1)/2$.
Using the Eistenstein integers, your equation can be factorized $$x^2 +3 y^2 = \Bigl(x+ (1+2 \omega) y \Bigr)\Bigl(x+ (1+2 \omega^2)y\Bigr)= N.$$
So given the unique prime-factorization of Eisenstein numbers (up to the units $\pm 1,\pm \omega,\pm \omega^2$), we have that $N$ should be factorizable in the same form.
So you can first factoring $N$ in the ordinary integers and thus we only need to discuss the case where $N$ is a prime number $p$. Now there are two cases
1) $p = -1 \pmod 3$: then $p$ is an Eisenstein prime and the equation has no solution (this is the case if any of the prime factors of $N$ are of this form).
2) in the other cases, we can factoring $p = (a + b \omega)(a+ b \omega^2) = a^2-ab +b^2$, we can then compare the prime factors to determine $x$ and $y$.
So the problem for the numerical algorithm is to solve $p=a^2-ab+b^2$ for $p$ a prime congruent to $1$ modulo 3 for which you know that a solution exists (the cases $p=3$ is covered by $a=1$, $b=2$). I have to admit that I do not know about any efficient algorithm for that except for trying $a$ and $b$.
Edit: here, is an algorithm which given a prime $p$ which is congruent to 1 modulo 3 finds $a$ and $b$ with $p= a^2 -a b +b^2$.
First, we note that $\gcd(a,b)=1$, that is $a$ and $b$ are coprime. All coprime integers can be generated by the Farey sequence where Wikipedia has an imlementation in Python. If the prime $p$ is not too large the algorithm of taking the coprime integers from the Farey sequence and simply testing if $p=a^2-ab+b^2$ will be efficient enough.