How to model interaction term with squared regressors? I run a regression with interaction of a squared continuous regressor with a categorical regressor. In Stata the double cross operator ## produces all combinations of my regressors. Here is an arbitrary MWE:
* load data
use http://www.stata-press.com/data/r13/nlswork

* set panel structure
xtset idcode year

* fixed effects regression with interaction and square term
quietly xtreg ln_wage c.wks_ue##c.wks_ue##i.race union age i.year, i(idcode) fe 
estimates store model1 

* generate square term without interaction
gen wks_ue_sq = wks_ue^2
quietly xtreg ln_wage c.wks_ue##i.race wks_ue_sq union age i.year, i(idcode) fe 
estimates store model2

estimates table model1 model2, keep(wks_ue c.wks_ue#c.wks_ue race#c.wks_ue race#c.wks_ue#c.wks_ue wks_ue_sq ) b p


Should I always go for what I call "full moderation"? As I would do in a threefold interaction. There is a theoretical reason why I have my square term, and for the interaction in general. But I cannot give a good argument about what the square-term-interaction would actually represents in my case.
 A: Two answers.  First, inclusion of the squared interaction is nearly always justified, in theory, by the Stone-Weierstrass theorem, which (loosely) states that your true mean response function is better and better approximated by higher order polynomials in your $X$ variables. The squared interaction term is one such higher order polynomial term. Of course, in practice, one must be concerned about overfitting, and about use of approximating polynomials because of their poor extrapolation properties.
However, this justification seems not so great in your case because you want a subject matter rationale. A second answer, albeit post hoc and therefore not ideal, is to estimate two models: (i) the simpler one, without the squared interaction but with all else deemed relevant, and (ii) the same model but including the squared interaction.  Then, construct profile plots of the estimated mean of $Y$ as a function of your continuous $X$ for fixed values of your categorical $X$. Overlay the profiles for model (i) and model (ii), given you one graph for each level of your categorical variable. Comparing these plots, along with use of your subject matter knowledge, should help you to decide what your squared interaction term is doing, and whether it is important.
Edit, 8/21/2020.  One thing that can be anticipated a priori when the squared continuous predictor interacts with the categorical predictor is that the character of the curvature will differ by categorical level. For example, there may be pronounced curvature within some levels but not others.
