What are methods I can use to calculate the two-sided Fisher exact test? I am implementing code to compute the two-sided Fisher exact test of 2x2 contingency tables in C. I read that R's implementation is not the most efficient way to calculate the p-value for a two-sided test. My first question is, how is this implemented? Wikipedia states
"An approach used by the fisher.test function in R is to compute the p-value
 by summing the probabilities for all tables with probabilities less than or
 equal to that of the observed table."

But I do not understand where the other tables and probabilities come from.
Also, how could I implement this function in a more computationally effective way, assuming that the statement
"A simple, somewhat better computational approach relies on a gamma function
 or log-gamma function, but methods for accurate computation of
 hypergeometric and binomial probabilities remains an active research area."

is correct?
 A: 
"An approach used by the fisher.test function in R is to compute the p-value
   by summing the probabilities for all tables with probabilities less than or
   equal to that of the observed table."
But I do not understand where the other tables and probabilities come from.

Fisher conditions on both margins. The other tables are all tables with the same margins. The probabilities come directly from the hypergeometric distribution. 
(Both these pieces of information are already in the article you point to, by the way.)
The help on the R function you mention mentions an algorithm for getting the other tables, but in the 2x2 case it's pretty simple, because they correspond to taking one cell in the table through the range of its possible values (the range of values taken by the random variable in a hypergeometric distribution). Here's one example (y-axis is log-scaled, but 
the values on the axis are actual probabilities):

The red numbers represent the ranks of the probabilities from smallest to largest. This 
is the order you account for them in computing a two-tailed p-value until you hit the observed table.

"A simple, somewhat better computational approach relies on a gamma function
   or log-gamma function, but methods for accurate computation of
   hypergeometric and binomial probabilities remains an active research area."

This is simply referring to the calculation of a hypergeometric probability. When the margins take large (or even middling) values, the exact calculation can be problematic (the calculations involve $\binom{n}{r}$ terms that can get quite large). It's usual to work on the log-scale with functions that return the log of the gamma function of its arguments.
R, for example, provides the functions lgamma and lfactorial for such purposes.
exp(lfactorial(5))
[1] 120
> exp(lgamma(6))
[1] 120

Fast methods for exact tests are often based on the approach of Mehta and Patel (variations of the network algorithm, mentioned, again at the link you gave, though there were more papers than listed there). There have been more recent developments in this area, but I am not especially familiar with them; if you know how to search for academic research (which these days mostly means "can google efficiently") more recent work won't be hard to find.
If you look at the references at the bottom of the help for R's fisher.test, it includes references to Mehta and Patel and also to the more recent paper by Clarkson, Fan and Joe which is more efficient. 
Page 17 of "Exact Inference for Categorical Data", 1997 by Mehta and Patel (on Cytel's website, here) mentions alternative approaches in passing (with references).
There's papers like 
Requena, F., & N.M. Ciudad, (2006)
"A major improvement to the Network Algorithm for Fisher’s Exact Test in 2 × c contingency tables,"
Computational Statistics & Data Analysis, 51, 490–498.
that might be of some use.
