How does one find the area of an implicit function? For example we have the equation $y^2+\sin({4y\cos{x}})=4$
You can see the graph here at:
https://www.desmos.com/calculator/1sxvfl2amd
So far I know it is split into top and bottom. I'm trying to find the area below the relation at the top, positive area region to the x-axis. Even if I isolate x to equal y or y to equal x, it still won't be a solvable integral.
I think I could compute the area numerically using the Riemann sum, but how can I use it specifically for the positive top, without isolating any variables. I'll fill in the limits from ${a<x<b}$, and/or $c<y<d$ of the integral to solve the rest myself, but I'm not sure how to?
 A: Since you are going to do the integral numerically, there is no loss in finding $y$-values numerically as well. Since the function you integrate is $2\pi$-periodic, and you want to integrate it over its period $[0,2\pi]$, there is no need for advanced methods like Simpson's rule or Gaussian quadrature; the simple left  endpoint rule will do just fine.
Here is for Sage, which can be used online (in a browser), free of charge.  It follows the logic of the Scilab code given below (which I commented in detail). 
var('n a b x y s') 
n = 1000
a = 0   
b = 2*pi
h = (b-a)/n
s = 0
for i in range(n):
    x = a + i*h
    s = s + find_root(y^2+sin(4*y*cos(x))-4, 1.5, 2.5)
float((b-a)/n*s)

The only difference, besides syntax, is that Sage asks for an interval in which to look for the root. I gave it $[ 1.5, 2.5]$, because your plot shows that the solution is in this range.

Here is code is for Scilab (requires installation, but is free), Text after // are my comments.
function z=f(y)
  z = y^2+sin(4*y*cos(x))-4    // implicit equation
endfunction

n = 1000      // number of sub intervals
a = 0         // left endpoint
b = 2*%pi     // right endpoint
h = (b-a)/n   // size of subinterval
p = a:h:b     // partition points, the last one (b) won't be used
v = []        // array to hold the values of function
for i = 1:n
    x = p(i)               // set x equal to partition point
    v = [v, fsolve(2,f)]   //  solve for y, taking y=2 as initial guess
end
integral = h*sum(v)              // left endpoint rule  
mfprintf(6, '%.9f', integral)    // formatted output, 9 digits

Output: $12.567055628$. I got the same answer for $n=2000$ subintervals, so it looks accurate.
Not coincidentally, this is close to $4\pi$, because the values of the function are close to $2$ on average. 
