Letting $x_1$ and $x_2$ be the planar coordinates, the most general cubic response has the form
$$\begin{aligned}
z & = \beta_0 \\
& + \sum_{i=1}^2 \beta_i x_i\\
& + \sum_{1\le i \le j \le 2} \beta_{ij}x_ix_j\\
& + \sum_{1\le i\le j \le k \le 2} \beta_{ijk}x_ix_jx_k.
\end{aligned}$$
The number of coefficients is $1 + 2 + 3 + 4 = 10.$ Although you propose $11$ locations, which ought to be enough to estimate all ten parameters, it's worth a check because one or more of those locations might be redundant for these purposes.
The way to check is to consider the design matrix for a single observation per location. If it is of full rank, the design matrix with multiple observations will be of full rank, too; and if it is not of full rank, replicating observations will not cure that.
Although we could do this check manually, it's so easy to get your statistical software to help that we might as well just put it to work. Begin with an array of coordinates of those 11 points:
df <- data.frame(x1 = c(0,0,0, 10,10, 20,20, 30,30, 40,40),
x2 = c(0,80,320, 40,160, 0,80, 40,160, 0,320))
Now pass it to your ordinary least squares routine, using any response vector you like. If it returns estimates for all $10$ variables, you're good. Here I used a vector of zeros for the response:
lm(rep(0, nrow(df)) ~ x1+x2 + x1*x2+I(x1^2)+I(x2^2) + I(x1^3)+I(x1^2*x2)+I(x1*x2^2)+I(x2^3), df)
Coefficients: (Intercept) x1 x2 I(x1^2) I(x2^2) I(x1^3) I(x1^2 * x2) I(x1 * x2^2) I(x2^3) x1:x2
0 0 0 0 0 0 0 0 0 0
Ten coefficients were estimated: your design is identifiable.
In R
, you can avoid computing and writing out all the terms by using its command to form polynomials, as in
X <- cbind(1, poly(as.matrix(df), 3, raw = TRUE))
This is a matrix with ten columns, one for each variable (including the intercept). A good way to analyze it is to compute its singular values:
print(svd(X)$d, digits=2)
[1] 4.7e+07 3.0e+06 7.8e+04 5.1e+04 1.7e+04 6.0e+02 2.2e+02 4.6e+01 2.9e+00 9.9e-01
They are all nonzero--the smallest, at the end, is 0.099--but there is a huge range. That's a sign you might have difficulties estimating some coefficients (or combinations) of those coefficients accurately.
R
(and many other platforms) offers a way to improve the numerical stability of the model through orthogonal polynomials. The correct analysis, as Belsley, Kuh, and Welsch demonstrate in Regression Diagnostics, is to analyze the centered variables (leaving out the intercept):
X <- cbind(1, poly(as.matrix(df), 3, raw = FALSE))
print(svd(scale(X)[,-1])$d, digits=2)
[1] 4.63 4.10 3.65 3.33 3.12 2.89 2.72 1.33 0.25
That's much better: the range is $4.63:0.25 = 18.7.$ This is the condition number of the design matrix. It's a little large: values exceeding $5$ in an experiment suggest room for improvement. As an example of what can be accomplished, this design has a much lower condition number
df <- data.frame(x1 = c(0,0,0, 10,10, 20,20, 30,30, 40,40),
x2 = c(0, 160, 320, 80, 240, 0, 320, 80, 240, 0, 320))

X <- poly(as.matrix(df), 3, raw = FALSE)
print(d <- svd(scale(X))$d, digits=2)
[1] 4.1 3.7 3.6 3.3 3.2 3.0 2.9 2.5 1.5
The singular values now run from $4.1$ down to $1.5$ for a condition number of just $2.8$. If this design meets your other criteria and is feasible to carry out, it will result in substantially more precise estimates of some coefficients for the same experimental cost (assuming this least squares model of the experimental variation is reasonable).
I am not proposing this particular redesign, but only pointing out how one can easily go about evaluating possibilities before going to the expense and trouble of doing the experiment.