How to solve the Gauss' law? Imagine we have the first Maxwell Equation:
$$
\nabla \cdot \mathbf E = \frac{\rho}{\varepsilon_0} \\
\mathbf E = -\nabla \phi \\
\nabla \cdot (-\nabla\phi) = \frac{\rho}{\varepsilon_0} \\\nabla^2 \phi = -\frac{\rho}{\varepsilon_0}
$$
We arrive to an equation for finding the potential of an Electric field, given a charge density, $\rho$.
Assuming that: $\phi$ and $\rho$ only depends on the $\mathbf r$, coordinate, we have the following:
$$
\nabla^2 \phi = \frac{\partial^2 \phi}{\partial r^2} + \frac{1}{r} \frac{\partial \phi}{\partial r}  = \frac{\rho(r)}{\varepsilon_0}
$$
For zero charge density we have the following:
$$
\frac{\partial \phi}{\partial r} \equiv u(r)\\
0 = u'(r) + \frac{u(r)}{r} \rightarrow u'(r) = - \frac{u(r)}{r} \\
\frac{du}{u} = -\frac{dr}{r} \rightarrow u(r) = \frac{C_1}{r} \\
\phi(r) = C_1\log(r) + \phi_0
$$
But, this mean there exists a potential. I don't quite understand this because there is not charge density, hence, no electric field. But how is there an electric potential, did I derive it wrong?
 A: You can have a non-zero potential -- as long as it is constant in space. This will always generate a zero electric field, since in this case $\mathbf{E} = -\nabla\phi = 0$, which is indeed expected for a zero charge density.
Mathematically, you could notice that since the point $r=0$ is part of your domain, you must have $C_1 = 0$ to avoid divergences in the potential. This generates $\phi = \phi_0$, i.e., a constant potential, as expected.
EDIT: As jensen paull said, you have the wrong Laplacian operator for spherical coordinates, making your final solution wrong. Nevertheless, as the correct solution is $\phi = C_1/r + \phi_0$, the above analysis still hold.
A: Ive been unfortunately trying to apply simple boundary conditions on your equation with no luck, because it's wrong. so First of all, let me correct you:
Laplaces equation for $r-$dependance is:
$$\frac{1}{r^2}\frac{\partial }{\partial r}\left(r^2\frac{\partial \phi}{\partial r}\right) = 0\implies r^2\frac{\partial \phi}{\partial r} = C_0 \implies \frac{\partial \phi}{\partial r} = \frac{C_{0}}{r^2}\implies \phi= \frac{C_{1}}{r} + \phi_{0}$$
Poissons equation for electrostatics.
$$\nabla^2 \phi = -\frac{\rho}{\epsilon_0}$$
There are 2 ways you can think about setting $\rho =0$ here.
Line of reasoning 1:
Lets say I have the equation $$z=yx$$
If I say $$y=1$$ This equation reduces to $$z=x$$
By substituting y=1, we are finding what z is, in regions where y=1.
That is not to say, that y=1 always. We are just finding the value of z that satisfy the condition that y=1
The same can be said for poissons equation.
$$\nabla^2 \phi = -\frac{\rho}{\epsilon_0}$$
Reduces to laplaces equation:
$$\nabla^2 \phi = 0$$
Using the former as a guide, what does this equation actually mean? By substituting $\rho = 0$ into poissons equation. We are saying "what equation does the potential $\phi$ follow, in regions where the charge density is zero?"
That is not to say that $\rho = 0$ everywhere, but where $\rho$ IS zero, we know the potential must follow this equation. This is true for all distributions. Take the laplacian of any potential in a region where there is no charge density, and the result is zero!
If poissons equation did automatically yield zero, then that would be even more confusing that if it didn't! As we KNOW the potential of e.g a point charge isn't zero, and most of the functions domain $\rho=0$.
Here the solution we find from laplaces equation ONLY works in regions where $\rho=0$, you may hear the phrase "finding the potential on a specific domain"
To find the constants we need to enter boundary conditions. Like the field is continuous about the boundary between the charged sphere/outside of it, or that the potential vanishes at infinity.
Line of reasoning 2:
We have discussed setting $\rho$ to zero to find the potential in regions where the charge density is zero. But what about if the charge density is actually zero everywhere?
The same equation should equally apply!
$$\phi= \frac{C_{1}}{r} + \phi_{0}$$
Given $$\phi(\infty) = 0$$
We then have $$\phi_{0} = 0$$
Given we know $\rho=0$ everywhere, we have no domain restriction, this fact will allow us to solve for $C_{1}$
We now need to say that $$\phi(0) ≠ \infty$$ Which is pleasing given we are saying there is no charge.
Doing so gives:
$$C_{1} = 0$$
Meaning $$\phi = 0$$
Non uniqueness of maxwells equations:
Maxwells equations are not uniquely determined by the curl and divergence of the fields.
Infact, the equation $\nabla^2 \phi = 0$ where $\rho$ is infact zero everywhere, does not need to be zero. The only reason that it is zero, is due to our boundary conditions, and assumption about the symmetry of our potential
In general, Any harmonic function can be added to the potential that does not need to be a constant, and will equally satisfy maxwells equations. This is the homogenous solution to maxwells equations in electrostatics. The helmholtz theorem states that if the fields vanish at infinity, then there is a unique solution. This is the normal  boundary condition that we impose on maxwells equations in order to find a unique solution, much like the field of a point charge.
Without it, there exists multiple   fields for the same charge distribution! For example the homogenous wave equation that can be added to any distribution, that is not uniquely determined by the charges, but by the boundary conditions of the EM field! [I wonder what the true boundary conditions are!!,pretty exciting that there could be a background EM field existing independantly of charges]
A: 
I don't quite understand this because there is not charge density, hence, no electric field.

This is wrong. There can be a non-zero electric field $\vec E(\vec r)$ at a point $\vec r$ even when the charge density at that point $\rho(\vec r)$is zero.
This is the case because the Laplacian $\nabla^2$ is a "non-local operator." This basically means just what I wrote in the sentence above. For example, a single isolated point change is a charge density that is only non-zero at a single point. But a single isolated point charge produces a field that is non-zero everywhere.

But how is there an electric potential,

There can be a non-zero potential because the potential is related to the charge by a non-local operator, as discussed above.

did I derive it wrong?

You did derive it wrong, but the wrongness does not really affect your question. You have the wrong numerical coefficient in one of your terms and you dropped a minus sign. You also used the $\partial$ symbol when the $d$ symbol can be used (since the $\phi$ only depends on $r$). The equation you wrote should be written as:
$$
-\frac{d^2 \phi(r)}{dr^2} - \frac{2}{r}\frac{d\phi(r)}{dr} = \rho(r)/\epsilon_0
$$
In the comments, you also asked:

How should I apply this conditions to have zero potential?

I take this to mean, how should you be able to tell when a charge density that is locally zero implies a potential that is also locally zero. The answer is that you have to solve the differential equation and therefore you must know the boundary conditions.
If the boundary conditions are such that $\phi=0$ on the entire boundary $\partial V$ of a volume $V$, then $\nabla^2\phi = 0$ (everywhere within the volume $V$) implies that $\phi=0$ everywhere within the volume $V$.
You can see this by considering the following positive-definite integral:
$$
\int_V d^3r |\vec\nabla \phi|^2 
= \int_V d^3r \left[\vec \nabla \cdot\left(\phi\vec\nabla\phi\right) - \phi\nabla^2\phi\right]\;.
$$
The second term on the RHS is zero (because $\nabla^2\phi=0$ everywhere in V).
The first term is also zero because:
$$
\int_V d^3r \nabla \cdot\left(\phi\vec\nabla\phi\right) = \int_{\partial V} d\vec S \cdot \left(\phi\vec\nabla\phi\right)\;,
$$
by the divergence theorem. But we also have that $\phi=0$ on the boundary $\partial V$, so the value of the integral is also zero.
Therefore:
$$
\int_V d^3r |\vec\nabla \phi|^2 = 0\;,
$$
which (since the integrand is positive definite everywhere) implies that
$$
\nabla \phi = 0
$$
everywhere in $V$.
So, at this point we see that $\phi$ has to be constant in $V$, but we can also see that the constant has to be zero because $\phi$ is zero at the boundaries.
You can also see this formally by considering the line integral from any point $\vec b$ on the boundary to any point $\vec r$ in the volume:
$$
0 = \int_{\vec b}^{\vec r} d\vec\ell \cdot \nabla\phi = \phi(\vec r) - \phi(\vec b) = \phi(\vec r)\;.
$$
