This question is related to one I asked over on physics stackexchange, but at this point it is purely mathematical in nature and I thought it would make sense to ask here. As the title suggests, I want to find an expression for the elements of a Jacobian, where the vector valued function in question is expressed as a Taylor series. More precisely I want to show that:
\begin{equation} \tag{1} \bigg(\frac{d \mathbf{P}^{(n)}}{d \mathbf{E}}\bigg)_{ij} = n \chi^{(n)}_{i jk_1 ... k_{n-1}} \prod_{i\in I_{n-2}} E_i \end{equation}
where $P^{(n)}$ is the n-th order term in: \begin{equation}\tag{2} P_i = \chi_{ij}^{(1)} E_i+ \chi_{ijk}^{(2)} E_i E_j+ \chi_{ijkl}^{(3)} E_i E_j E_k + ... \end{equation}
and $I_n$ is a set containing the indices of the n-th order tensor. The tensor $\chi^{(n)}$ has $n+1$ Cartesian indices, and is symmetric with respect to permutation of these indices. I know how to show this for say the 2nd order term:
\begin{equation}\tag{3} \bigg(\frac{d \mathbf{P}^{(2)}}{d \mathbf{E}}\bigg)_{ij} = 2\delta_{jk} \delta_{jm} \chi^{(2)}_{imk}E_k+ (1-\delta_{jm}) \delta_{jk} \chi^{(2)}_{imk}E_l + (1-\delta_{jk}) \delta_{im} \chi^{(2)}_{imk}E_k \end{equation}
In the first term, because $m=k$, there's a factor 2 when taking the derivative. In the second term I can simply relabel the indices, and then use the permutation symmetry to cancel all the terms with two delta functions, leaving:
\begin{equation}\tag{4} 2\delta_{jm}\chi_{imk}^{(2)} E_k = 2\chi_{ijk}^{(2)} E_k \end{equation}
Similarly for the 3rd order term:
\begin{equation}\tag{5} \begin{gathered} \bigg(\frac{d \mathbf{P}^{(3)}}{d \mathbf{E}}\bigg)_{ij} = 3\delta_{jm}\delta_{jk}\delta_{jl} \chi_{imkl}^{(3)}+ 2(1-\delta_{jm})\delta_{jk}\delta_{jl} \chi_{imkl}^{(3)} E_m E_k\\+ 2(1-\delta_{jk})\delta_{jm}\delta_{jl} \chi_{imkl}^{(3)} E_m E_k+ 2(1-\delta_{jl})\delta_{jk}\delta_{jm} \chi_{imkl}^{(3)} E_l E_k\\+ (1-\delta_{jm})(1-\delta_{jk})\delta_{jl} \chi_{imkl}^{(3)} E_k+ (1-\delta_{jl})(1-\delta_{jm})\delta_{jk} \chi_{imkl}^{(3)} E_k\\+ (1-\delta_{jl})(1-\delta_{jk})\delta_{jm} \chi_{imkl}^{(3)} E_l \end{gathered} \end{equation}
All terms with more than one delta-function cancel, leaving:
\begin{equation}\tag{6} 3\delta_{jm}\chi_{imkl}^{(2)} E_kE_l = 3\chi_{ijkl}^{(2)} E_kE_l \end{equation}
from the last three terms in the long sum above. Now, I don't know how to show this for the general case, or even a good way how to express it. I could define something like:
$$\tag{7} \delta_{ij}^{\mu} = \begin{cases} \delta_{ij} & \mu = 1 \\ (1-\delta_{ij}) & \mu = 0 \end{cases} $$ and then I should take the product of $n-1$ such functions:
\begin{equation}\tag{8} \prod_{j\in I_{n-1}} \delta_{ij}^{\mu_{j}} \end{equation}
but summed over every combination of $\mu_{j}$, so something like:
\begin{equation}\tag{9} \sum_{\mu_1=0}^{1}...\sum_{\mu_{n-1}=0}^{1} \prod_{j\in I_{n-1}} \delta_{ij}^{\mu_{j}} = \prod_{j=0}^{n-1} \sum_{\mu_j = 0}^{1} \prod_{j\in I_{n-2}} \delta_{ij}^{\mu_{j}} \end{equation}
Each product of delta functions should also be multiplied by the sum of $\mu_j$ to account for the factor from the derivative:
\begin{equation}\tag{10} \prod_{j=0}^{n-1} \sum_{\mu_j = 0}^{1} \sum_{j=1}^{n-1}\mu_{j}\prod_{j\in I_{n-1}} \delta_{ij}^{\mu_{j}} \end{equation} In total:
\begin{equation}\tag{11} \bigg(\frac{d \mathbf{P}^{(n)}}{d \mathbf{E}}\bigg)_{ij} = \prod_{j=0}^{n-1} \sum_{\mu_j = 0}^{1} \sum_{j=1}^{n-1}\mu_{j}\prod_{j\in I_{n-1}} \delta_{ij}^{\mu_{j}} \chi^{(n)}_{i_1 ... i_n} \prod_{i\in I_{n-2}} E_i \end{equation}
but I'm not convinced this is correct, and I also don't know how to usefully manipulate that expression at all. I'm sure there's some better way of approaching this problem that I'm not seeing.
Nonlinear optics! One of my favorite subjects.
Anyway consider the following gradient:
$$\frac{\partial (\prod_{i\in I_n} E_i)}{\partial E_k} = \sum_{j\in I_n} \frac{\partial E_j}{\partial E_k} \prod_{i\neq j}{E_i} = \sum_{j\in I_n} \delta_{jk} \prod_{i\neq j} E_i = n\prod_{i\in I_{n-1}} E_i $$
where the last equality holds by the symmetry of the arguments (permutation of the indices). Given that $\chi^{(n)}$ shares that same symmetry (and is constant w.r.t. $\mathbf{E}$), I believe the relation you want follows.