Why does electric field intensity $E$ can be uniquely determined by its divergence and curl? My question is, the number of following equations
$$\nabla\cdot E=\frac{\rho}{\varepsilon}$$
$$\nabla\times E=-\frac{\partial B}{\partial t}$$
is 4 while the number of unknown variables $E=(E_1,E_2,E_3)$ is 3. Intuitively, the equation is overdetermined and the solution may not exist unless four equations are correlated. Is my intuition right?
 A: The divergence and curl do not uniquely determine a vector field. For example, if all the derivatives are zero, any constant field is a solution. In order to determine the field, we must include boundary conditions.
As far as your question, you are not looking at equations for the three components of $E$. You are looking at equations for the partial derivatives of those components, of which there are nine ($\frac{\partial E_x}{\partial x}, \frac{\partial E_x}{\partial y}$), etc.
A: The electric field $E$ is uniquely determined by its curl (Faraday's law) and its divergence (Poisson's) if it is assumed that the field must fall to zero at infinity.  The two operators are indeed correlated so that they do not overdetermine the field.
To see what's going on, suppose $E_2$ is a solution of the second (curl) equation (Faraday's law).  Then so is $E_1+E_2$, where $E_1=-\nabla  \phi$ is a gradient of a scalar function.  Any scalar function (within reason). (In math parlance, gradients are in the kernel of the curl operator.)  
So, you need 1 more equation to determine the scalar function $\phi$.  And, voila, there it is, the first equation (Poisson's).  
To be tidy about it, move all the gradients (stuff with non-zero divergences) into the $E_1$ field, so that $E$ is the sum of two components with one ($E_1$, known as the longitudinal component) curl-free and the second ($E_2$, aka the transverse component) divergence-less.  Qmechanic provides explicit formulas for these components in his answer to this question.
Now, given such a field $E$ has been constructed, are there other fields which also satisfy the two equations (i.e. have the same divergence and curl)?  By the Helmholtz decomposition theorem, this solution is unique, if it is required that the field falls off to zero at infinity "sufficiently quickly".  (You can see the precise definition of "sufficiently quickly" in the reference, as well as the required reasonableness of the functions involved.)  For example, you can add a constant field and the equations are still satisfied.
All this assumes that:
1) $B$ is known and is divergence-less
2) the currents $j$ are such that the last Maxwell equation is satisfied.

Update.  More on uniqueness:  If $E_a$ and $E_b$ have identical divergences and curls, their difference $E_d$ has zero curl (and is hence the gradient of a scalar function $\psi$ in a simply connected domain) and zero divergence, so $E_d$ is the gradient of a harmonic function (solution of Laplace's equation): 
$$E_d=-\nabla \psi \quad , \quad  \nabla^2 \psi = 0$$ 
Assuming $E$ goes to zero at infinity forces $E_d=0$ (and proves uniqueness), via a Green's Theorem application.  See for example here (pdf).
A: Let's look at this from a numerical perspective.  Take Maxwell's equations and solve them for the time derivatives.  (Here, I use $\epsilon_0 = \mu_0 = c = 1$ units.)
$$\begin{align}
\nabla \cdot E &= \rho \\ 
\frac{\partial E}{\partial t} &= \nabla \times B - j \\
\nabla \cdot B &= 0 \\
\frac{\partial B}{\partial t} &= - \nabla \times E 
\end{align}$$
From initial data on some spatial hypersurface, the time evolution of the EM field is determined by the 6 time derivative equations (2) and (4) above.  The other two equations involve only spatial derivatives and form constraints on the initial data.
For general differential equations of this form, I think the takeaway here is that boundary conditions are constrained.  Let me restrict myself to the case $\nabla \cdot E = \rho$ and $\nabla \times E = 0$.  We know that the electric field at any point can be found using a convolution with the Green's function--that is, the field of a point charge--along with a surface integral:
$$E(r) = \int_M \frac{\rho(r') (r-r')}{4\pi|r-r'|^3} \, dV' + \text{nasty surface integral over $\partial M$}$$
Usually the surface integral is ignored because we take the boundary condition at infinity, but it is definitely there, and it involves integrating $E$ directly on that surface (along with the Green's function).
But once $\rho$ is chosen, $E$ is constrained on the surface to respect the choice of $\rho$.  It must because we have, through Gauss's law,
$$\int_M \rho dV = \oint_{\partial M} E \cdot dS = Q$$
In other words, you cannot simply put any choice of $E$ on the surface as a boundary condition.  You can put some $E$ that respects the total charge within the volume, and then perhaps add any extra $E'$ that corresponds to charge from outside the volume (and thus is divergence- and curl-free), but that's all the freedom you have.
In short: boundary condition data is constrained.  The additional constraint equations eliminate degrees of freedom and should make everything match up.
A: Maxwell's equations do not uniquely determine $\mathbf{E}$ and $\mathbf{B}$.
Instead, they place constraints on how the fields may coexist given different setups. Usually, one must provide boundary conditions to get the actual answer; in physical situations there are almost always boundary conditions which bring  back the deterministic nature of the system.
This is rather evident when you attempt to solve the equations for a moving rod in a time varying (but space-constant) magnetic field. Quoting from my answer to the relevant question

Your fundamental issue is that Maxwell's equations (of which Faraday's
  law is one) are not "cause and effect". You cannot "plug in" a value
  of magnetic field and get a corresponding value of $\bf E$ field
  induced by the $\bf B$ field. All Maxwell's equations tell you is
  "which kinds of $\bf E$ and $\bf B$ fields can coexist given so-and-so
  conditions".
Trying to solve the situation via Maxwell's equations
I remember solving a similar situation via Maxwell's equations and
  being surprised by the answer.
The "initial conditions" were $\mathbf {B}=\beta t\hat k$, $\rho=0$
  (no charge), $\mathbf{J}=0$ (no current).
Solving{*} for $\mathbf{E}$, using the differential+microscopic form
  of Maxwell's equations(since the integral form can only get you the
  value of $\bf E$ at certain positions at many times), I got:
$$\mathbf{E}=\hat i (lx + \frac{\beta}{2}y+az+c_1)+\hat
> j(-\frac\beta{2}x+my+bz+c_2)+\hat k(ax+by+nz+c_3)$$
where $a,b,l,m,n,c_1,c_2,c_3$ are arbitrary constants subject to
  $l+m+n=0$
Note that this is a family of electric fields (Setting certain
  constants to zero, you get concentric ellipses IIRC). All this means
  is that any $\bf E$ field of this type can coexist with a $\bf B$
  field.
Implication for your problem
This means that your initial conditions are insufficient/inconsistent.
  Along with such a magnetic field, any type of electric field
  satisfying the above equations can exist--and must exist. 
So, in addition to knowing how your magnetic field is changing with
  time, you need to know:
  
  
*
  
*Which one of these bajillion electric fields is present
  
*Where is the rod in relation to this electric field?
  

