@shoyer Attached the description of the issue source and kind of workaround.

During concat a result_vars = OrderedDict() is created. After that it is iterated over the first dataset datasets[0].variables.items() and those variables which are not affected by the concat are added to the result_vars : https://github.com/pydata/xarray/blob/a5ca64ac5988f0c9c9c6b741a5de16e81b90cad5/xarray/core/combine.py#L244-L246

After several checks the affected variables are treated and added to result_vars: https://github.com/pydata/xarray/blob/a5ca64ac5988f0c9c9c6b741a5de16e81b90cad5/xarray/core/combine.py#L301-L306

The comment indicates what you already mentioned, that the reorder might be unintentional. But due to the handling in two separate iterations over datasets[0].variables, the source variable order is not preserved (and with that in some cases the order of the dimensions).

This can be worked around by changing the second iteration to: ```python

re-initialize result_vars to write in correct order

result_vars = OrderedDict()

stack up each variable to fill-out the dataset (in order)

for k in datasets[0].variables: if k in concat_over: vars = ensure_common_dims([ds.variables[k] for ds in datasets]) combined = concat_vars(vars, dim, positions) insert_result_variable(k, combined) else: insert_result_variable(k, datasets[0].variables[k])
`` With this workaround applied, theconcat` works as expected and the variable/coordinate order (and with that the dimension order) is preserved. I'm thinking about a better solution but wanted to get some feedback from you first, if I#m on the right track. Thanks!