In the $(k, d)$-Function-SAT problem we are given a set of $n$ variables $\{X_1, \dots, X_n\}$ that can take values from the set $\{1, \dots, d\}$ and a set of Boolean constraints on these variables, where each constraint is of the form $f:\{1, \dots, d\}^k \right arrow \{0, 1\}$, i.e. the constraint depends on exactly $k$ of these variables. We will treat $k$ and $d$ as constants. The goal is to determine whether the set of constraints has a satisfying assignment, i.e. an assignment to the variables such that all constraints simultanuously map to $1$. In this paper, we study $(k, d)$-Function-SAT in the property testing model for dense instances. We call an instance $\epsilon$-far from satisfiable, if every assignment violates more than $\epsilon n^k$ constraints. A property testing algorithm is a randomized algorithm that, given oracle access to the set of constraints, must accept with probability at least $3/4$ all satisfiable inputs and rejects with probability at least $3/4$ all inputs, which are $\epsilon$-far from satisfiable. We analyze the canonical non-adaptive property testing algorithm with one-sided error: Sample $r$ variables and accept, if and only if the induced set of constraints has a satisfying assignment. The value of $r$ will be called the \emph{sample complexity} of the algorithm. We show that there is an $r_0 = \wide tilde O(1/\epsilon)$ such that for any instance that is $\epsilon$-far from satisfiable, the probability, that a random sample on $r\ge r_0$ variables is satisfiable, is at most $1/4$. This implies that the above algorithm is a property tester. The obtained sample complexity is nearly optimal for canonical testers as a lower bound of $\Omega(1/\epsilon)$ on the sample complexity is known. Previously, a tester with sample complexity $o(1/\epsilon^2)$ was only known for the very special case of testing bipartite ness in the dense graph model \cite{AK02}. Our new general result improves the best previous result for testing satisfiability (and even for the special case of $3$-colorability in graphs) from sample complexity $\wide tilde O(1/\epsilon^2)$ to $\wide tilde O(1/\epsilon)$. It also slightly improves the sample complexity for the special case of bipartite ness. Improving the sample complexity for $(k, d)$-Function-SAT (or special cases of it) had been posed in several papers as an open problem \cite{AK02, AS03, RS11}. This paper solves this problem nearly optimally for canonical testers and, in the case of $k=2$, also for non-adaptive testers as there is a lower bound of $\Omega(1/\epsilon^2)$ on the query complexity of non-adaptive testers for bipartite ness in the dense graph model \cite{BT04}, where the query complexity denotes the number of queries asked about the graph (for a canonical tester in graphs, the query complexity is the square of its sample complexity). As a byproduct, we obtain an algorithm, which, given a satisfiable set of constraints, computes in time $O(n/\epsilon^{O(1)} + 2^{\wide tilde O(1/\epsilon)})$ a solution, which violates at most $\epsilon n^k$ constraints.