I feel a bit out of my condifent zone here, but playing around, I think the following is rough but gives some insight. taking the previous code for granted.
What if we define a separate varaible z, to capture the higher order relationship of the poly x order 2.
If we fit the gam with x +s(z) where s(z) is some smoothed x^2 pure function with no linear term...
I can compare the plot(gam2) which is the y~x+s(x) with my hand plot of x by an approximate f(x,order 2) function) and I think its a close enough match.
Hopefully someone else can way in.
d$z <-((d$x-300)^2)/150 -200
(gam3 <- mgcv::gam(y ~ x + s(z),data=d,method="REML"))
(s3 <- summary(gam3))
par(mfcol=c(2,2))
plot(gam3,all.terms = TRUE)
par(mfcol=c(1,2))
plot(gam2)
plot(d$x,d$z)