Many　physical　processes　appear　to　exhibit　fractional　order　behavior　that　may　vary　with　time　or　space.　The　continuum　of　order　in　the　fractional　calculus　allows　the　order　of　the　fractional　operator　to　be　considered　as　a　variable.　Numerical　methods　and　analysis　of　stability　and　convergence　of　numerical　scheme　for　the　variable　fractional　order　partial　differential　equations　are　quite　limited　and　difficult　to　derive.　This　motivates　us　to　develop　efficient　numerical　methods　as　well　as　stability　and　convergence　of　the　implicit　numerical　methods　for　the　space-time　variable　fractional　order　diffusion　equation　on　a　finite　domain.　It　is　worth　mentioning　that　here　we　use　the　Coimbra-definition　variable　time　fractional　derivative　which　is　more　efficient　from　the　numerical　standpoint　and　is　preferable　for　modeling　dynamical　systems.　An　implicit　Euler　approximation　is　proposed　and　then　the　stability　and　convergence　of　the　numerical　scheme　are　investigated.　Finally,　numerical　examples　are　provided　to　show　that　the　implicit　Euler　approximation　is　computationally　efficient.